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Preliminary results concerning the time evolution of strongly exited SU{2) Bogomolny-Prasad- 
Sommerfield (BPS) magnetic monopoles have been published in [13] . The behavior of these dynam- 
ical magnetic monopoles was investigated by means of numerical simulations in the four dimensional 
Minkowski spacetime. The developed code incorporates both the techniques of conformal compact- 
ification and that of the hyperboloidal initial value problem. Our primary aim here is to provide 
a detailed account on the methods and results of the investigations reported in In addition, 

some important new results, which go much beyond the scope of these early studies, are also pre- 
sented. In particular, to be able to distinguish linear and non-linear effects, evolutions of monopoles 
deformed by various excitations, including both very small and extra large energy excitation, are 
investigated. In addition, a detailed account is provided on the spacetime dependence of the basic 
variables, as well as, that of the other physically significant dynamical qualities such as the energy 
and energy current densities, the radial and angular pressures and the magnetic charge density. A 
careful comprehensive study of the associated energy transfers and energy balances is also included. 

PACS numbers: 03.50.Kk, 14.80. Hv 



INTRODUCTION 



Soliton and quasi-soliton type configurations, which 
are in general spatially localized non-singular finite 
energy solutions of non- linear field equations, play 
significant role in various particle physics considera- 
tions (see e.g. jH, [131 for recent reviews). In general 
the associated field equations are highly non-linear, 
therefore, most of the investigations have been re- 
stricted to the study of time independent configura- 
tions. Thereby, there is an obvious increase of in- 
terest to study the dynamical properties of soliton 
type configurations. It is a matter of fact that the 
available analytic techniques either have not been de- 
veloped enough to provide or simply cannot guaran- 
tee completely satisfactory answers to all questions in 
case of non-linear systems. Therefore, reliable numer- 
ical approaches are needed which are able to describe 
the time evolution of these type of systems. Moti- 
vated by this sort of necessities we have developed 
a numerical method which have been used in study- 
ing the time evolution of various spherically s vmmet 
ric non-linear dynamical systems [ill, [H, O, [ 
also 0, [^, [^, 01 )• The numerical approach we ap 



(see 



plied is based on the "method of lines" in a fourth 
order setting jl^l which was found to be the most 
efficient among the numerical methods which are ap- 
plied in various considerations [l^l- In addition, to 
get a faithful representation of all the radiation pro- 
cesses, the techniques of conformal compactification 
(introduced first by Penrose [3l|) to the underlying 
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Minkowski spacetime, along with the hyperboloidal 
initial value problem were incorporated. The applied 
conformal gauge, which is a modification of the static 
hyperboloidal gauge, makes it possible to investigate 
by numerical means the asymptotic properties of ra- 
diation processes. One of the associated advantages is 
that by making use of this method the radiation pro- 
cesses, which last infinitely long in physical time, can 
be analyzed in finite computational time intervals. 

Among the great variety of physically interest- 
ing soliton type configurations distinguished attention 
has been paid to the study of magnetic monopole type 
configurations such as the 't Hooft-Polyakov mag- 
netic monopole solutions of coupled Yang-Mills-Higgs 
(YMH) systems [32J. However, even in the case of 
magnetic monopoles there was a lack of knowledge 
concerning the dynamical properties of these systems. 
Therefore, our numerical method was first applied 
to the study of the dynamical properties of the sim- 
plest possible magnetic monopole system, which could 
be chosen so that there was a radiative component 
among the basic field variables. This paper is, in fact, 
to report about the results of investigations concern- 
ing the time evolution of a strongly excited spherically 
symmetric SU{2) BPS magnetic monopole Pf. 

To have a clear enough setting, at least from an 
"energetic" point of view, we investigated the time 
evolution of an initially static SU(2) BPS magnetic 
monopole which was excited by the help of a high en- 
ergy pulse. The self-interaction of the Higgs field was 
turned off which yielded a Yang-Mills-Higgs system 
so that the Yang-Mills field is massive while the Higgs 
field is massless. 

The dynamics started by injecting a pulse of excita- 
tion via the time derivative of the Yang-Mills variable. 
Since the Yang-Mills and Higgs variables are coupled 
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(in a non-linear manner) energy was transformed to 
the massless Higgs field which immediately forwarded 
about half of the energy of the pulse towards future 
null infinity, J^'^ . It is already a surprise that this 
reaction of the monopole is apparently independent 
of the energy content of the exciting pulse, at least 
in case of the subclass considered in this paper (see 
subsection I VII B 2]) . Nevertheless, naively, one might 
still expect that the following part of the evolution is 
going to be quite boring since only a lower scale en- 
ergy transform will happen until the rest of the energy 
of the exciting pulse disperse and the system settles 
down to the static monopole. 

Contrary to this simple hypothetical scenario quite 
interesting features of the underlying non-Hnear sys- 
tem show up. First of all, it was found that the ex- 
citing pulse leaves the central region so that it drags 
some part of the energy of the original static monopole 
with itself. This process is justified by the fact that 
in the central region, where a long-lasting quasi-stable 
'breathing state' develops, at least at the beginning of 
the evolution the average energy density is always less 
than that of the static monopole was. Since the orig- 
inal static magnetic monopole is stable we know that 
by the end the missing energy has to come back to the 
central re gion . However, as our numerical simulations 
indicates this process is unexpectedly slow and, 
in certain cases, behaves completely contrary to our 
expectations. For instance, the difference between the 
dynamical and the original static energies contained 
in a ball of radius r may be positive or negative de- 
pending on the strength of the non-linear aspects of 
the dynamics, as well as, on location, i.e. on the value 
of r. In particular, the time dependence of this "ex- 
tra energy" content can be characterized by its high 
frequency oscillating part and by its mean value. We 
found a "universal" time decay for the amplitude of 
the high frequency oscillations of the form of a power 
law with exponent close to the "ideal" value —5/6. 
However, the exponent of the similar power law time 
decay of the mean value varies with the radius of the 
ball (see subsection I VII C 5|) . 

It is also a characteristic feature of the evolution 
that some part of the energy of the exciting pulse is 
getting to be stored by expanding shells of oscillations 
of the massive Yang-Mills field built up in the distant 
region. These shells of oscillations behaves much like 
the oscillations of a simple massive Klein-Gordon field 
(see [HI for a detailed investigation) since the cou- 
pling of the Yang-Mills and Higgs variables is negligi- 
ble there. The monopole can pull back the missing en- 
ergy, approaching the static configuration in the cen- 
tral region, only before these shells get too far away. 

Due to the coupHng of the Yang-Mills and Higgs 
fields in the central region, along with the "breath- 
ing" of the monopole there, the massless Higgs field 
continuously take away a small fragment of the en- 
ergy of the monopole, which is radiated to future null 



infinity. The corresponding decrease of the total en- 
ergy which can be associated with the hyperboloidal 
hypersurfaces decreases in time with power —2/3 (see 
subsection lVII C 4]) . Note, however, that this total en- 
ergy (it would be called as the 'Bondi energy' in gen- 
eral relativity), as well as, its limit value, is always 
strictly larger than the energy of the original static 
monopole. This is because the expanding shells of os- 
cillations store some part of the energy of the original 
pulse forever. In fact, this part of the energy cannot 
reach future null infinity since it is associated with 
the massive Yang-Mills field, moreover, as they are 
getting further and further away from the central re- 
gion the monopole can get back less and less energy 
from them since non-linear effects are getting to be 
more and more neghgible. As it is shown in subsec- 
tion IVII C 5| the energy stored by these expanding 
shells also approaches its non-zero asymptotic value 
from above following a power law time decay with 
exponent —2/3. 

The time evolution of the frequency of oscilla- 
tions of the excited monopole behaves also somewhat 
counter to the general expectations. Consider first the 
fundamental frequency of the "breathing" monopole. 
Instead of having a slow decay in the frequency of the 
associated oscillations as their amplitude decreases 
there is an increase in the frequency which takes the 
maximum value at the 'end of time', i.e. at future 
timelike infinity. Interestingly, while the frequency 
of the breathing monopole tends to the mass of the 
Yang-Mills field from below the frequency of the os- 
cillations stored in the expanding shell has no up- 
per bound, although the amplitude is decreasing, as 
in case of the massive Klein-Gordon field Note 
that these phenomena are not at all unphysical. Re- 
call that a simple physical system Hke the spinning 
coin on the desk can produce a quite similar effect. 
As energy loss happens due to dissipative processes 
the amplitude of the "spinning" is getting smaller and 
smaller in consequence of which the frequency of the 
oscillations is getting higher and higher. In case of 
the distant oscillating shells the unbounded increase 
of the frequency is, in fact, to compensate the de- 
crease of the oscillation amplitudes which together 
ensure the conservation of the energy stored by the 
shells. 

In virtue of these results it seems not to be over- 
rating to say that the investigation of the dynamics 
properties of magnetic monopoles made transparent a 
number of interesting and unexpected features of the 
underlying non-linear system. Hopefully, these sort 
of investigations will stimulate further numerical and 
analytic investigations of various similar non-linear 
dynamical systems. For instance, partly motivated 
by the findings of the above described investigations, 
in a framework of linear perturbation theory some 
of the above mentioned features could be explained 
successfully both qualitatively and quantitatively 
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Note, however, that truly non-Hnear effects are by 
their very nature out of the scope of these sort of 
investigations. In particular, as some of the new 
results of this paper make it transparent, they cannot 
be used to describe the behavior of the investigated 
system when the energy of the exciting pulse is much 
larger than that of the original static monopole. 

Although in this paper we use our numerical code to 
describe various oscillations of a localized monopole, 
we would like to note that oscillations in a non- 
linear theory may be responsible for the very ex- 
istence of some compact objects. In a wide vari- 
ety of nonlinear field theories, where no static soli- 
ton solutions exist, almost periodic long living oscil- 
lon configurations have been found. See for example 
EHHHIIll,!!!! for recent results on this quickly 
evolving topic. A slightly modified version of our nu- 
merical code have also been used in for a detailed 
study of oscillon configurations in cj)'^ scalar theory. 

The structure of this paper is as follows. In the 
next section, after recalling some of the basics related 
to the properties of the underlying generic dynamical 
system a detailed description of the specific choice for 
both the YMH system and the geometry is provided. 
The static hyperboloidal conformal gauge applied in 
our investigations is introduced in section llllj while 
the first order representation of the field equations, 
relevant for the used conformal setting, is given in 
section IIVI A detailed description of the appHed nu- 
merical scheme is presented in section [Vj while vari- 
ous numerical tests of the code, for the case of massive 
and massless Klein-Gordon field, are presented in sec- 
tion IVII Finally, all the numerical results concerning 
the evolution of dynamical magnetic monopoles, in- 
cluding subsections providing detailed description of 
the time dependence of the basic and derived field 
variables and that of the frequency of the associated 
oscillations, the energy transfers, monitoring of the 
numerical violation of constraints and the energy con- 
servation, along with the behavior of the radial and 
angular pressures, as well as, the magnetic charge den- 
sity, are presented in section IVIll 



II. PRELIMINARIES 

The investigated dynamical magnetic monopole is 
described as a coupled SU{2) YMH system. The 
Yang-Mills field is represented by an su(2)-valued vec- 
tor potential Aa and the associated 2-form field Fab 
reads as 



Fab = daAb - dbAa + ig [Aa, Ab] 



(1) 



where [ , ] denotes the product in su(2) and g stands 
for the gauge coupling constant. The Higgs field (in 
the adjoint representation) is given by an su(2)-valued 
function ip while its gauge covariant derivative reads 



as Vail = daip + igiAajip]- The dynamics of the in- 
vestigated YMH system is determined by the action 

J {Tr{F,fF-f) + 2 [Tr(I?eV'2?'^) " e, 

(2) 

where e is the 4-dimensional volume element, more- 
over, describing the self-interaction of the Higgs 
field, is chosen to be the standard quadratic potential 

V{^) = ^[Tr{i?)-Hl]\ (3) 

where A and Hq denote the Higgs self-coupling con- 
stant and the 'vacuum expectation value' of the Higgs 
field, respectively. 

The symmetric energy-momentum tensor of the 
considered YMH system takes the form 



Tab — —-r- 



TriFaeFb") - Tr{Va^Vb^) + ]^gabC 



(4) 



where C stands for the Lagrangian 

C = Tr{FefF''f) + 2 [Tr(X>eV'2?''V') - Vii')] ■ (5) 

A. Fixing the geometrical and gauge setup 

This paper is to investigate the evolution of spher- 
ically symmetric Yang-Mills-Higgs (YMH) systems 
on a fiat Minkowski background spacetime. Accord- 
ingly, as a fixed background, the four-dimensional 
Minkowski spacetime (R^,?yah) will be applied, the 
line element of which in the conventionally used 
Descartes-type coordinates (a;° ^) reads as 

ds^ = {dx°f - {dx^f - {dx^f ~ {dx'^f. (6) 

The gauge group is specified by giving the set of 
generators {r^-} (1=1,2,3) of the associated su(2) Lie 
algebra which reads as 



1 



(7) 



where cTj denote the Pauli matrices. The commuta- 
tion relations relevant for this choice of generators are 



(8) 



where e^,^ denotes the completely antisymmetric 
tensor with e^^^ — 1. 

In addition to the above special choice concerning 
the gauge group, our considerations will be restricted 
to YMH systems which are yielded by the 'mini- 
mal' dynamical generalization of the static t'Hooft- 
Polyakov ma gne tic monopole configurations [13, lUl 
(see also (l^. l25l|). Accordingly, the Yang-Mills and 
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Higgs field variables, Aa — A'^t^ and — tp'T^^ are 
specified, in the Coulomb gauge, via the relations 



ij' = H- 



(9) 



(10) 



where the functions w and H are assumed to be 
smooth functions in depending upon the coordi- 
nates , x^,x^,x^ only in the combinations t — x'^ 
and r = \/{x^y + (x'^Y + {x^y. 

Since the metric, as well as the matter fields, are 
required to be spherically symmetric, the use of the 
standard coordinates {t, r, 9, (/)), adapted to the spher- 
ical symmetry of {W^,r]ab), is the most suitable. In 
these coordinates the line element of the Minkowski 
metric takes the form 



+ sin^ 9 d,p'^). (11) 

We have assumed above that w and H are smooth 
functions of t and r. This might be surprising espe- 
cially because our eventual aim is to carry out numer- 
ical simulations of the YMH systems under consider- 
ations. Numerical methods are inherently too rough 
to make a sensible distinction between configurations 
belonging to different differentiability classes. In fact, 
our smoothness assumption is to ensure certain tech- 
nical conveniences used later and it is supported by 
the following considerations. (Note that for the fol- 
lowing argument it would be enough to assume that 
w and H are smooth in a sufficiently small neighbor- 
hood of the origin.) In numerical simulation of spher- 
ically symmetric configurations a grid boundary rep- 
resenting the origin inevitably appears. This means 
that we need to solve an initial-boundary value prob- 
lem which, in particular, requires the specification of 
'boundary behavior' of w and H at the origin through- 
out the time evolution. This is, in fact, the very point 
where we make use of our smoothness assumption. In 
virtue of ^ and (flOl) it is straightforward to see that 
spherical symmetry, along with the required smooth- 
ness of w and H through the origin, ensures that in 
a neighborhood of the origin w and H are even and 
odd functions of the r-coordinate, respectively. Con- 
sequently, as it is described in subsection lV Dl in more 
detail, by extending our grid by a suitable number 
of virtual grid points the 'boundary values', for in- 
stance, of the r-derivatives of w and H need not be 
specified by hand, instead they are naturally yielded 
by the time evolution of the corresponding odd and 
even functions. 

The Yang-Mills field Aa as it is given by ^ is in the 
"Coulomb gauge". Another frequently used gauge rep- 
resentation is the so called "abelian gauge" in which tp 
has only one non- vanishing component, and which is 
achieved by making use of the gauge transformation 



A straightforward calculation yields that in the corre- 
sponding abelian gauge the Yang-Mills and the Higgs 
fields read as 

A..-i[.{r.(d.)„-r,sin.(d0)4 + r3COS.(d^).] 

(13) 

and 



(14) 



The substitution of these gauge representations into 
the equations of motion, deducible from the action 
|[2|) , yields the evolution equations for w and H which 
are given as 



r^d^w - r^dtw = w [{w^ - l) + g^r^H^] (15) 



r^d^H + 2rdrH - r^d?H 



H 



2w^ + -^r^iH^ - Hi) 



(16) 



The static finite energy solutions of these equations 
are called 't Hooft-Polyakov magnetic monopoles. 



B. Regularity at the origin and at infinity 

These equations, along with the former assumption 
concerning the smoothness of the basic variables, en- 
sure certain regularity of the Yang-Mills and Higgs 
fields. In particular, it turned out that the values of 
w and H are restricted at the origin and at spacelike 
infinity throughout the evolution. At the origin (see 
appendix A for more details) the relations 



w{t,Q) = 1, and H{t,Q) = 



(17) 



have to hold to have finite energy density, measured 
by the static observer with four velocity u°- = (d/dt)"". 
In addition, whenever the smoothness of w and H is 
guaranteed in a neighborhood of the origin we also 
have 



di:w\r^o = 
for k being odd, while 



(18) 



(19) 



U = cxp(z0r2 ) exp(z(; 



(12) 



for k being even. For the first sight it might be a bit 
of surprise but a same type of reasoning that lead to 
lfT8|) and lfT9|) provides a restriction also for the rest 
of the derivatives. Namely, from the field equations 
and from the smoothness requirement it follows that 
at the origin all the non-vanishing derivatives d!pw, 
with k > 4, and dlH, with I > 3, can be given as 
a function of d^w and drH, along with their various 
time derivatives. 
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The argument associated with the regularity of the 
basic variables at spacelike infinity is more delicate. 
To get it we need to refer to the conformal (non- 
physical) setting, where our numerical simulation ac- 
tually will be carried out. In this framework it seems 
to be essential to assume that the fields are at least 
even through future null infinity, J^'^ . We refer 
here to this requirement as being merely a technical 
assumption, nevertheless, we would like to emphasize 
that the use of it is supported by the following obser- 
vations. In [s^l it was proved that the time evolution 
of a massive Klein-Gordon field in Minkowski space- 
time with initial data of compact support necessar- 
ily yields C(l/r°°) asymptotic behavior at null infin- 
ity. In addition, the results of Eardley and Moncrief 
0, [3, [13, [131 concerning the local and global existence 
of YMH fields in 4-dimensional Minkowski spacetime 
supports that the above technical assumption can be 
deduced from the field equations at least in case of 
suitably chosen initial data specifications. 

This smoothness requirement, along with the rele- 
vant form of the rescaled field equations, implies (see 
appendix A for more details) that w has to vanish 
while H has to be constant along the null geodesic 
generators of . The field equations also yield fur- 
ther restriction on this constant limit value of H 
at J'^ . In particular, must take the value 
whenever A 7^ but it is an arbitrary (positive) con- 
stant, C, otherwise. Thus, by making use of a limiting 
argument, the values of w and H at spacelike infinity 
can be seen to be determined as 



lim w{t, r) — 0, 



(20) 



H = lim H{t,r)^l^'' ^ ^ ^' (21) 

^ r-»oo I otherwise. ' 

Notice that the later relations are in accordance with 
our general expectations that the considered dynam- 
ical YMH systems do really possess the asymptotic 
fall-off properties of a magnetic monopole [la . [25| 
throughout the time evolution. 

It follows from equation ifTSj) . lfT6|) and (|2T|) that 
the semi-classical vector boson and Higgs mass, 
and M„ , are 



(22) 
(23) 



In particular, both fields are massless whenever 
vanishes, while only the Higgs field is massless, al- 
though it is coupled to a massive vector boson, when- 
ever the self coupling constant A is zero. 

Whenever the smoothness of the basic variables 
is guaranteed through then the n*^ order r- 
derivatives of w and H are also restricted there. 
In particular, it follows from the relations ijA.lOp . 
(|A.15P - l|A.18p . along with the vanishing of the energy- 
momentum expressions at null infinity that whenever 



the asymptotic fall off conditions, 

lim \r^''d^w{t, r)] = (for n > 0) (24) 

r — >oo 

and 

lim [r^" 9;' iJ(i,r)] = 0, (for n > 2) (25) 

— ^00 

have to be satisfied. 

It is important to note that whenever A 0, i.e. 
the self-coupling of the Higgs field is turned down, 
there exists an explicitly known [s^l static solution to 
ifTSl) and iflGl) 



gCr 



= C 



sinh(gCr) 
/ 1 



1 



\tanh((7Cr) gCr 



(26) 
(27) 



where C is an arbitrary positive constant which is 
in fact the limit value, H^, of the r.h.s. of (|27| . 
This static solution is called to be the Bogomolny- 
Prasad-Sommerfield (BPS) magnetic monopole which 
is known to be a linearly stable configuration [T|. 

Notice that in the BPS limit, i.e. whenever A = 
and ~ Hq ^ 0, the Higgs field becomes mass- 
less and the only scale parameter of the system is the 
vector boson mass — gH^ . Since in the case con- 
sidered here 7^ 0, the re-scahngs t i = tM^, 
r ^ f = rM^ and H ^ H = H/H^ transform the 
parameters to the value g — = 1. This implies 
that whenever we would like to study the excitations 
of the BPS monopole it suffices to consider the time 
evolution of the system only for the particular choice 
of the parameters g = = 1, because by inverse 
re-scalings of these solutions all the possible solutions 
to the field equations can be generated. Actually, to 
check the efficiency of our numerical implementation 
of the evolution equations first the stability of this 
static solution was investigated. Later, the complete 
non-Hnear evolution of systems yielded by strong im- 
pulse type excitations of these analytic static solution 
was also studied in detail. 

Although it is not necessary to be done, neverthe- 
less we preferred to put the principal part of (fT6|) into 
the same form as the principal part of ifTSj) which was 
achieved by making use of the substitution 



H{t,r) 



h{t, r) 



H 



(28) 



It follows then that our new basic variable h vanishes 
at the origin while it takes a finite limit value at infin- 
ity if H is guaranteed to tend to its boundary value 
fast enough there. The substitution of (|28|) into the 
relations (fT5|) and iflG]) yields then 



- 9» = 
{h 



w WW 



l)+g\h + H^ 



2w^ 



A 



\{h^ + 2H^rh) 



(29) 



,(30) 
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where we have used the relation = Hq when- 
ever A ^ 0. The regularity of the solutions to these 
equations at r = follows from the boundary con- 
ditions (flTl) . It is important to note that l(2T|) and 
(|28l) . along with our smoothness and symmetry as- 
sumptions, which in particular implies that H is an 
odd function of r, guarantees that h + H^r has to be 
an even function of the r-coordinate. Moreover, since 
h? + 2H^rh = {h + H^rf - H'^^r'^ the right hand 



sides of l(29|) and ((30|) are both even functions of the 
r-coordinate. Notice finally that, by equations (|29l) 
and l(30l). the expressions d'jw — dfw and d^h — dfh 
have regular limits at r = oo because, in virtue of 
(EU, 1211), {111 and ([281, w and h remain finite while r 
tends to infinity. 



III. THE STATIC HYPERBOLOIDAL 
CONFORMAL GAUGE 

In this section a conformal transformation of the 
Minkowski spacetime - along with the relevant form 
of the above matter field equations - will be consid- 
ered. 

The applied conformal transformation is a slightly 
modified version of the gauge transformation applied 
first by Moncrief (s^l (see also [IB, [13). It is defined 
by introducing the new coordinates (T, R) instead of 
(i, r) as 



T{t, r) — Kt — \J K^r^ + 1 
1 



R{v) 



1 



with the inverse relations 



l-i?2 



and 



2R 



K(l-i?2)' 



(31) 

(32) 
(33) 



where n is an arbitrary positive constant. The 
Minkowski spacetime is covered by the coordinate do- 
main satisfying the inequalities — oo < T < +oo and 
< i? < 1 (see FigEl). 

The line element of the conformally rescaled metric 
Tjab = ^^gab in coordinates (T, R, 9, (f>) takes the form 



^2 



2RdTdR - dR^ 

- i?2 [dO"^ +siiv^ed(j}'^) , (34) 



where the conformal factor is 

K 



n{R) = -{i-R'), 



(35) 



and, by ((33|) and ([35|) . we have the relation 

rVl = R. (36) 



In this conformal representation the R = \ coor- 
dinate line represents J^'^ through which the metric 
Ijab smoothly extends to the coordinate domain with 
R>1. 

The name 'static hyperboloidal gauge' is explained 
by the following observations. First, l(32|) tells us that 
the R — const lines represent world-lines of 'static ob- 
servers', i.e. integral curves of the vector field {d/dt)°-. 
Second, it follows from (|3T|) that the T = const hy- 
persurfaces are, in fact, hyperboloids satisfying the 
relation {Kt — T)^ — n^r'^ = 1 in the Minkowski space- 
time. 

The equation describing radial null geodesies in the 
(T, R) coordinate system is independent of the param- 
eter K, 



T : 



±2R 
l±R 



To 



(37) 



with the plus signs corresponding to outgoing and the 
minus signs to ingoing geodesies. The form of the 
constant is chosen in a way that both type of geodesies 
cross the origin i? = at time T = Tq. The outgoing 



Outgoing radial null geodesies --.^ \ \\\z 

"Ingoing" radial null geodesies \ lijs 

















0.6 
R 



FIG. 1: The ingoing and outgoing radial null geodesic curves 
are shown in both the "physical" and "unphysical" regions. 
As in other conformal representations each point represents a 
two-sphere of radius R. The interior of the strip < i? < 1 
corresponds to the original Minkowski spacetime while J^"*" is 
represented by the R = 1 coordinate line. Notice that the R = 
const lines possesses timelike character everywhere except at 
the R = I line which is null. 

geodesic emanating from the origin reaches future null 
infinity i? = 1 in a finite coordinate time at T = Tq-I-I. 
The ingoing geodesies starting from a point at radius 
i?o at T = reaches the origin at T = 2i?o/(l — 
i?o)- Consider now a geodesic starting at a point close 
to J^"*" with radius Rq — 1 — 6R. Such a geodesic 
reaches the origin at T = — 2 -I- 2/SR. This means 
that for geodesies coming in from the far away region, 
the travel time essentially doubles when the geodesic 
starts from a point "twice as close" to the R = 1 
line. This has an important consequence in relation 
to the validity domain of our numerical code. Namely, 
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the doubling of the resolution yields the doubling of 
the time interval within which we may expect our 
numerical simulation to provide a proper solution to 
the selected problem. Here it is assumed tacitly that 
all the possible inaccuracies born at and coming from 
the outer region - that is the part of the spacetime 
where the evolution cannot be described properly by 
any of the numerical techniques based on a finite grid 
- travel inwards only with the speed of light. 

Consider now the field equations relevant for this 
conformal setting. To start off take a function / = 
f{t,r) of the coordinates t and r and denote by 
/ = /(T,i?) the function f{t(T,R),r{T,R)) yielded 
by the substitution of (|33| into /. A straightforward 
calculation justifies then that for any pair of functions 
/ = f{t, r) and / — f{T, R) the relation 



-2RdRdTf 



+ 1)2 \^ ^2 R 



(38) 



29. 



K(i?2 + 1) 



holds. It follows then from ^ and ^ that the 
field equations, l(29|) and l(30|) . in the current confor- 
mal representation read as 



(i?2 + 1)2 

4i?2 

(i?2 + 1)2 



'^W — w 



l)+g^(h + H^Rn- 



1 \2 



h + H RQ 
2w^ + - 



h{h + 2H RQ- 



(39) 



,(40) 



where the second order partial differential operator *P 
is defined as 



02 

_/92 
«;2'^« 



- 2RdRdT 



20 



K(i?2 



QRiR^ + 3) 
K(i?2 + 1) 



(41) 



As before, the regularity of the solutions to these 
equations at the origin, i.e. at i? = is guaranteed 
by the regularity of |[29]) and |[30]) at r = 0, while, the 
regularity at future null infinity, i.e. at i? = 1, follows 
from |[2l]), ([Ml) and ([Ml) provided that the fields have 
at least extensions through J^^ . 

Hereafter, unless indicated otherwise, all the func- 
tions will be assumed to depend only on the coordi- 
nates T and R. Thereby we suppress all of the tildes, 
introduced above. 



IV. THE FIRST ORDER HYPERBOLIC 
SYSTEMS 

This section is to derive a first order hyperbolic 
system from the above evolution equations for which 



the initial value problem is well-posed. To see that 
equations l(39|) and l(40|) can be put into the form of 
a strongly hyperbolic system we shall follow a stan- 
dard process (see e.g. [^). Correspondingly first we 
introduce the first order derivatives of w and h, 

wt = drw, wr = Orw, hx = drh, Hr = dRh, 

(42) 

as new variables. In terms of the relevant enlarged 
set of dependent variables (|39)) and l(40l) can be given 
as 



drh 



02 

k2 

O^ 

2 



drWT = — {Orwr) -2R{dRWT) 

{dRhR)^2R{dRhT)' 



K 



- K 



(43) 
(44) 



where bj^ and bh are given as 

20 Oi?(i?2 



bw = 



K(i?2 + 1) 

(i?2 + 1)2 

w 

402i?2 

20 

~ k{R^ + 1) 

(i?2 



Wt 



3), 



-Wr 



k{R^ + 1) 
02 (^2 _ l)+ g2 (^fifi + JJ^Jl-) 

nR{R^ + 3) , 



(45) 



1) 



402i?2 



K(i?2 + 1) 

{nh + H^R) 

X 



(46) 



20^2 + -h{nh + 2H^R) 



It is straightforward to see that the above equa- 
tions, along with the first and third equations of (|42l) . 
and the integrability conditions 



Qtwr — Orivt and drhR = OrHt, 
possess the form of a first order system 

- AR{dR<P) + B, 
where $ and B are given as 



moreover, Ar takes the form 



(47) 
(48) 



/ w \ 




/ WT \ 


Wt 




bw 


Wr 

h 


, B = 




Ht 






bh 


KhRj 




\ J 



(49) 



/O 

-2R 02/^2 

1 





Vo 



-2R 02/k2 
1 



(50) 



Since the eigenvectors of Ar comprise a complete sys- 
tem and its eigenvalues are all real this first order 
system is, in fact, a strongly hyperbolic system [2^ . 
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Notice also that B does not depend on the space- 
hke derivatives of $, i.e. it is a functional of T, R and 
<i> exclusively, B = B{T,R;^). It is also straight- 
forward to verify that the constraint equations, i.e. 
the second and the fourth equations of l|42p . are pre- 
served by the time evolution governed by (|48l) . This, 
in particular, guarantees that the pair of functions w 
and h yielded by the evolution of suitable initial data 
specifications will automatically satisfy the original 
evolution equations l(39|) and l(40|) as well, provided 
that the constraint equations hold on the initial data 
surface. 

Finally, we would like to emphasize that the sys- 
tem specified by the relations l(48|) - (|5Q|) is not only 
strongly hyperbolic but can be put into the form of 
first order symmetric hyperbolic system. To see this 
introduce the variables wr and hn by the relations 

n - n 

Wr = —wu and hj^ = —hu- (51) 

K K 
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FIG. 2: The domains above the initial data hypersurface 
T = To covered by the simplest possible orthogonal grid are 
indicated in the applied conformal representations. 



Then the system of field equations for the vector vari- 
able $ = {w, WT, WR, h, hr, fiR)^ , possess the form of 
(I48l) with 



/O \ 

-2i? n/K 

n/n 



-2R n/n 

\o n/K / 



(52) 



Nevertheless, in all of our numerical simulations the 
strongly hyperbolic form given by l|48p - (|50|) was used 
which system was found to be as efficient in accuracy 
in a number of particular cases as the symmetric hy- 
perbolic system. 



boundaries at the origin and beyond future null infin- 
ity (see subsection lV DI for more retails). For this pur- 
pose we always used 1/2^ and 1/2'* parts of the spatial 
gridpoints on the left side of the origin and beyond 
future null infinity. This means that the formal val- 
ues of the i?-coordinate at these timelike portions of 
the numerical grid boundaries are Ri — 



-0.017 and Rr 



1- 



16 



256-4-16 

1.068. Note that the 



256-4-16 

appropriate treatment of the grid boundaries does not 
require the use of fixed portions of the spatial grid- 
points. Nevertheless, we used this simple approach 
to have a straightforward setting in which we could 
compare the results of our numerical simulations for 
various resolutions for the price that only the 236/256 
part of the spatial gridpoints represented points from 
the original Minkowski spacetime. 



THE NUMERICAL SCHEME 



In our numerical simulations we shall use the sim- 
plest possible orthogonal grid based on the T = const 
and R — const 'lines' in the domain T > To as it is 
indicated on Figl2l The relevant discrete set of grid 
points is given as 

Ti - To + IAT, l^O,l,...,L^ax (53) 

R, = iAR, i = 0, 1, (54) 

for some Imax , Lmax S N fixed numbers and with the 
relation AT — kAR for some k G N. In this setting a 
function / = /(T, R) will be represented by its values 
at the indicated grid points, i.e. by = /(T;, Ri). 

In our numerical simulations the total number of 
spatial grid points was chosen to be an integer power 
of 2, i.e. 2^, for some value N gN. A small fragment 
of these spatial gridpoints were used to handle the grid 



A. The time integrator 

The time integration of (|48l) is based on the use of 
the 'method of fines' in a higher order scheme as it is 
proposed by Gustafsson et al [23|- In particular, we 
integrate l(48|) . along the constant Ri lines, by mak- 
ing use of a fourth order Runge-Kutta scheme. This 
is done exactly in the manner as the Runge-Kutta 
scheme is used to integrate first order ordinary dif- 
ferential equations. Correspondingly, the value of $ 
after a 'time step' is determined as 
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where 



^ AT 



= AT 



(56) 



^B{Ti + -AT,Rr,<P\ + -^') , (57) 



= AT 



AT 



Ri', 

+B{Ti + AT,Rr,^[ + 'i',,,,,[)], (59) 



i3(r, + iAr,i?,;$', + ivl,^„,;)],(58) 



To be able to apply this method of time integration 
we must calculate i?-derivatives of certain functions 
several times. In this time integration process we ap- 
proximated these i?-derivatives by a symmetric fourth 
order stencil (see appendix B). 

It follows from general considerations ji^l that the 
above described time integration process applied to 
(|48l) cannot be stable unless a suitable dissipative 
term is added to each evaluation of the right hand 
side of l(48|) . An appropriate dissipative term, rele- 
vant for the fourth order Runge-Kutta scheme used 
here, reads as 



D =a(Ai?)5(a^$), 



(60) 



where cr is a non-negative constant and the sixth or- 
der i?-derivatives were evaluated in a symmetric sixth 
order stencil (see appendix B). 

Notice that the use of this dissipative term does 
not reduce the order of accuracy of the applied finite 
difference approximation. Moreover, since l(48|) is an 
almost everywhere symmetrizable strongly hyperbolic 
system and all the differential operators are centered 
it follows from Theorem 6.7.2. of f23| that the applied 
time integration process is stable provided a is suf- 
ficiently large and k = AT/AR is sufficiently small. 
Numerical experiments showed us that, for instance, 
the particular choices a ~ 10~^ and fc = 1/8 provide 
the required stability for the time integration of our 
evolution equation. 



B. The evaluation of the evolving fields at the 
origin and at future null infinity 

Before being able to make the first time step, and 
later at each of the time levels we have to face the fol- 
lowing technical problem. The source term B in (|48l) 
contains terms - see the last terms of l(45l) and (|46|) - 
which are given as ratios of expressions vanishing at 
the origin and at conformal infinity, respectively. The 
evaluation of these type of expressions generally is a 
hard numerical problem. 



It turned out, however, that the exact value of 5„ 
and bh is either irrelevant or can be shown to be zero. 
To see this consider first the evolution of the fields 
at the origin. By making use of the field equations 
and assuming that the fields are at least of class 
at i? = we have shown (see (fTTl - ^T9\\ ) that in 
a sufficiently small neighborhood of i? = the field 
values w and h must possess the form 



w 



1 + ?«(r, R) ■ R^ 



h = hiT,R)-R^-^H 



(61) 
(62) 



where w{T,R) and h{T,R) are assumed to be suffi- 
ciently regular functions of their indicated variables. 
These relations immediately imply that 



moreover, 



u;(r,0) = l, 
M;T(r,0) = 0, 
WRiT,0) = 0, 



hiT,0) = 0, 

/iT(r,o) = o, 



hn{T,0) = - 



2H 



(63) 
(64) 
(65) 



(66) 
(67) 

(68) 



hold. In other words, the time evolution of our basic 
variables at the origin is trivial and hence the evalu- 
ation of the bw and bh is not needed there. 

Despite their apparent singular behavior the evalu- 
ation of bw and 6;i at i? = 1 is possible based on the 
following observations. The first two terms of 5„ and 
bh are proportional to the conformal factor O so they 
vanish at . The last terms also vanish there in 
spite of the presence of the il"^ factors in them since 
these expressions always contain as a multiplying fac- 
tor at least one of the massive field variables, either w 
or, whenever A ^ 0, /i. These field variables, however, 
possess the fall off property that r^w and r"/i for arbi- 
trary positive integer value (n G N) tend to zero while 
r ^ Qo. By making use of this fact, along with the 
relation r = R/VL, it is straightforward to check that 
both bw and bh must vanish at J^^ . Accordingly, the 
values of 6^ and bh were kept to be identically zero 
at i? = 1 in the numerical simulations. 



C. Increasing accuracy 

In order to ensure higher order of accuracy we also 
have applied the following trick. Instead of calculat- 
ing the time evolution of the variables w and h them- 
selves we determined the evolution of the deviation, 
= w — wq and = h — hg, of them with respect 
to certain analytic functions, wo and ho. In fact, the 
functions wq and Hq need not necessarily required to 
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be solutions of the field equations although it is fa- 
vorable to assume that they possess the same type 
of behavior at i? = and at as the functions w 
and h themselves. This way it was possible to achieve 
a considerable decrease of the error of our numerical 
scheme in these critical neighborhoods. 

The evolution equations for and can be de- 
duced immediately by making use of the assumptions 
that w = -I- Wo and h — + ho satisfy i|48p . The 
linearity of the involved differential operators yields 
that 



The boundary at the origin 



where 



dT^^^AR{dR<i>^)+B^, (69) 



B^=B- 9t$o + AR{dR<i>o) (70) 



with B being the functional as it is given by l|49p but 
evaluated at the functions w = + wq and h = 
+ ho, moreover, <i>^ and $o denote the vectors 



(«'a)t 



/ ^0 \ 

dTWQ 

OrWo 
ho 
drK 
V dnho j 



(71) 



It is straightforward to verify that in the particular 
case when wo and ho are solutions of l|39p and l(40|) . 
i.e. $0 is a solution of l(48|) . 



B. 



B-Bo 



( (w'a)t \ 


bh - bho 




(72) 



where bw and bh are supposed to be evaluated at the 
functions w = + wo and h = h^ + ho, respec- 
tively. In all our numerical simulations, presented in 
this paper, wo and ho were chosen to be the static 
BPS solution given by ^ and ^ with C = g = 1. 



D. The treatment of the grid boundaries 

In any numerical simulation the appropriate treat- 
ment of the grid boundaries requires the most care- 
ful considerations. In the case studied here the grid 
boundary consists of two disjoint parts, the gridpoints 
representing the origin and the right edge lying in the 
non-physical spacetime region. It is the evaluation 
of various derivatives of the dependent variables at 
these boundaries what desire very careful investiga- 
tions. The handling of these issues, as they are im- 
plemented in our code, are described in the following 
subsections. 



As it was already mentioned in Section Bl Al at the 
origin we can use spherical symmetry and the smooth- 
ness of the field variables to gain information about 
the parity of the relevant functions. This knowledge 
can then be used to extend these functions onto a 
virtually enlarged grid. In particular, for any func- 
tion which can be divided into two parts, one of 
which has a definite parity in coordinate R, while 
the other part is explicitly known, suitable symme- 
try transformations can be used to define values of 
the function at the additional grid points to the left 
of the word line representing the origin (see Fig. [3]). 
By making use of these values the required number 






R 







FIG. 3: The way of extension of functions to i? < is 
indicated for even, and odd, functions, respec- 

tively. In fact, the extensions are determined by either an 
axial or a central reflection to the -R = line or to the 'central 
point' as indicated by the formulas {jeuen jea;t|_^ _ jeuen 

and {f°'^'^y^\-i = respectively. 

of i?-derivatives of the applied variables can be deter- 
mined at R — 0. More explicitly, as we will see any 
function / in our interest can be given as a sum of two 
functions, / — " + /" , one of which has a defi- 
nite parity /'' '' , being either even or odd function of 
the i?-coordinate, and another explicitly known one 
^ . Such a function can be extended, its extension 
will be denoted by {f}'^'^*, to the additional gridpoints 
as follows. Suppose, first, that f'^ p- is an even func- 
tion. Then f"^"^ = / — can be extended to the 
left of i? = simply by taking the axial refiection, 
|je„e„}ext ^ (/ „ f'-)^^, of its graph to the R = 0- 
axis, defined by the relation (/="'" )"^^|-i = f^^"". 
Finally, by making use the fact that /■= is explic- 
itly known and thereby its extension is supposed to 
be also known (potential singularities at the origin 
has already been excluded by former assumptions), 
the extension {f}'^^* of / can be given as {/}^^' — 
if ~ f"^ )^^- + A similar process appHes 

whenever f"^ " is an odd function with the distinction 
that axial refiection should be replaced by a 'central 
refiection' determined by the relation (/'"''' )*^^|-i = 
_jodd Thereby the extension of / then can be given 

as {/r^* = )^« + {/--r"*. 



11 



Turning back to our concrete field variables note 
first that in the time integration process we need 
to evaluate i?-derivatives of the dependent variables 
listed as components of the vector valued functions <&, 
see the relation l(49|) . Since w and h + rH^ were found 
to be even functions of the original r-coordinate, it 
follows by a straightforward substitution that 



w 
Wt 

h 



w 



W 
h 



R 



2172 



H 



(73) 
(74) 
(75) 

(76) 

(77) 

(78) 



Thereby the extensions of these functions, to the en- 
larged grid, yielded by the above described general 
process are 







(79) 






(80) 




= -wbI, 


(81) 




^M.2H.(5);, 


(82) 






(83) 


flB-i 




(84) 



Note that whenever we consider the evolution of 
A$ with respect to a solution $0 of l(48l) - $0 pos- 
sesses the same parity properties as the hypothetical 
solution $ itself - thereby the components of A$ are 
necessarily even functions of the i?-coordinate so their 
extensions at the origin are straightforward. 



2. The grid boundary in the non-physical spacetime 

The appropriate treatment of the boundary in the 
non-physical spacetime requires completely different 
type of considerations. Clearly, at this part of the 
boundary there is no way to enlarge the grid based 
on a suitable combination of certain smoothness and 
symmetry requirements, as it was possible to be done 
at the origin. Instead to be able to determine the re- 
quired number of i?-derivatives we used the following 
two ideas: 

First, at the edge and neighboring gridpoints we ap- 
plied a numerical adaptation what would be called as 
'one sided derivatives' (see appendix B) in analytic in- 
vestigations. By making use of these approximations 
the time integration process can be carried out along 
the gridlines close to the edge. Second, our problem 



inherently is a boundary initial value problem (see e.g. 
[2^). Thus appropriate boundary conditions have to 
be chosen at the grid boundary to be able to carry 
on the time integration scheme. There is a consider- 
able freedom in choosing boundary conditions at the 
edge of the grid. We intended to choose the one which 
allows waves to travel from the left toward the edge 
of the grid without being reflected, moreover, which 
excludes waves coming from beyond this edge toward 
the direction of the domain of computation. The cor- 
responding restrictions in terms of our basic variables 
in the static hyperboloidal gauge are 



Wt 



-^{R+lfwB., 



hT = --{R+lfhB. 



(85) 
(86) 



Numerical experiments justified that this choice did 
really ensure that the above two requirements were 
satisfied. 



E. Specific ciioices of parameters in the code 

In order to concentrate sufficient number of grid 
points to the central region where the monopole lives, 
as well as, to have enough grid points to resolve the 
expanding shell structures at large radius we have 
chosen k = 0.05 as the parameter included in the 
coordinate transformation ((311) and l(32|) . This spe- 
cific choice turned out to be appropriate for all the 
simulations presented in this paper. 

For convenience, the total number of spatial grid 
points was always chosen to be an integer power of 2. 
Our minimal resolution was 2^ — 256, with 4 points in 
the negative R region and 16 points in the R> 1 do- 
main. To preserve the grid point positions when dou- 
bhng the resolution, in general, we used 4/256 = 1/64 
part of the total spatial grid points for the negative 
R mirror image points, and 16/256 = 1/16 part of 
the total points for the unphysical domain. Keeping 
the "size" of the unphysical region above i? > 1 to 
be constant is useful when investigating the stability 
and the convergence of the code in that region, while 
the unnecessary points in the negative R domain only 
decrease the speed of the code slightly. 

In order to investigate the stability conditions of the 
numerical code it is helpful to write out the coordinate 
velocity of the radial null geodesies. Using l(37|) we get 



^ = ±i(l±i?)2 
dT 2^ ' 



(87) 



where the positive signs correspond to outgoing, the 
negative signs to ingoing geodesies. The maximal co- 
ordinate velocity occurs at null infinity i? = 1, tak- 
ing the value 2 for outgoing geodesies. The absolute 
value of the coordinate velocity of ingoing geodesies 
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is always smaller than 1/2. In principle this would 
allow any numerical time step AT which is smaller 
than half of the radial step AR. However, the evolu- 
tion equations of the massive field components have 
an apparent singularity at i? = 1, where our numeri- 
cal evolution code was not appropriately stable unless 
the time step was chosen as small as AT — AR/8. 

Choosing the parameter a in the dissipative term 
((60l) to be cr = 0.01 stabilize our numerical code with- 
out influencing the results signiflcantly even at lower 
resolutions. 



VI. TESTING THE NUMERICAL CODE 
WITH MASSIVE AND MASSLESS 
KLEIN-GORDON FIELDS 

In order to ascertain the appropriateness of our nu- 
merical code we employed it to a physical system for 
which the time evolution can be established by an in- 
dependent, certainly more precise method. Because 
of the linearity of the equations describing the mas- 
sive Klein-Gordon fleld its time evolution can be cal- 
culated by the Green function method. If the Klein- 
Gordon field and its derivative are given on a spacelike 
hypersurface then the field value can be calculated at 
any point in the future as a sum of two definite (nu- 
merical) integrals (for more details see, e.g., [12l])- The 
study of the Klein-Gordon field is especially impor- 
tant, since for several physically important systems, 
including the magnetic monopole, at large radius the 
various field components decouple, and behave like 
independent massive or massless linear Klein-Gordon 
fields. 

In the spherically symmetric case the Klein-Gordon 
field $ satisfies the wave equation 



Introducing the new field variable z = r$, in the co- 
ordinate system (T, R) defined by ^ and ^ the 
Klein-Gordon equation takes the form 



{R^ + 1)- 



m z 



(89) 



where the differential operator *P is defined in l(4T|) . 
Similarly to the monopole case the equation can be 
transformed to a system of three first order partial dif- 
ferential equations by introducing the new dependent 
variables zt = drz and zr = Brz. 

The initial data on the T — hypersurface was 
chosen as a specification of nonzero time derivative in 
a localized region superimposed on the vacuum value 
of the field z = as 



{zt)o 



.exp 



, if r G [a — 6, a + 6] ; 



0, 



otherwise, 

(90) 



with the constants selected to be a = 2, 6 =1.5, 
c — 70, c? = 10 and k = 0.05 in this section. 

First, we applied our evolution code to the massive 
Klein-Gordon field with m — 1. On Fig.|4] we com- 
pare the field values on a constant T slice obtained 
with different spatial resolutions to the precise value 
calculated by the Green function method. The chosen 
time slice is at T = 2.4746, well after the null geodesic 
emanating from the origin at T = has reached null 
infinity at T = 1. It is apparent, that even if the 

r 

150 200 300 400 600 1000 2000 





600 1000 2000 




FIG. 4: The behavior of the massive Klein-Gordon field 
with m — 1 aX the T = 2.4746 slice for large radiuses is 
shown. In order to make the expanding shell structures 
more apparent the function 2: = r<l> is plotted instead of 
On the upper graph the oscillations of z as calculated 
by the Green function method is shown. On the lower 
graph the upper envelope of the oscillations in z, i.e. the 
curve connecting the maximum points of the function z, 
is shown for the Green function calculation, and also for 
the evolution code corresponding to the indicated spatial 
resolutions. 

results obtained by the evolution code are not com- 
pletely satisfactory close to (i? = 1) the field val- 
ues in the central area remain correct for considerably 
longer time intervals. On Fig. [5] the upper envelope 
curve for $ is shown at a fixed radius, R — 0.02542, 
corresponding to r = 1.0176. The frequency of the 
field oscillations approaches (from above) m = 1 in 
terms of the original time coordinate which in the 
T coordinate corresponds to 20. We can see, that as 
it may be expected from the behavior of the ingo- 
ing null geodesies, the doubling of the resolution in- 
creases the time interval where the numerical solution 
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50 100 150 200 250 300 350 400 450 500 
T 

FIG. 5: The upper envelope curve of the oscillations of 
the Klein-Gordon field $ at a constant radius close to the 
center of symmetry. In order to compensate for the fast 
decay of the field, the envelope of $T^/^ is plotted, which 
should tend to a constant value on theoretical grounds. 
Apart from a short initial period, the envelope of the 
Green function result is really constant. The envelopes 
of the functions obtained by the evolution code stay near 
this constant value for longer and longer times as the 
number of spatial grid points is increasing. 



is valid approximately by a factor of two. On Fig. [6] 
the absolute value of the error at the same constant 
R radius is shown for several different numerical res- 
olutions. Initially, up to approximately T = 0.1, the 
error values decrease according to the expected fourth 
order convergence of the code. Later, for T > 0.3 the 
convergence become even faster, approximately fifth 
order. This behavior is due to dominance of the arti- 
ficial dissipation term which decreases as (AR)^ . At 
later stages a numerical instability arises and the error 
increases until its absolute value reaches the magni- 
tude of (j) and the numerical simulation is no longer 
reliable. 

We have also tested our code by applying it to 
the massless Klein-Gordon field with the same pa- 
rameters in the initial data ((90l) . With m = the 
Klein-Gordon equation ((89l) is regular at null infinity 
i? = 1 in the compactified representation. On Fig. [7] 
we show the rescaled field variable z = r$ at null 
infinity as a function of the hyperboloidal time co- 
ordinate T obtained with several different numerical 
resolutions. Since in the massless case information 
spreads strictly with the velocity of light, the value 
of z has to be exactly zero before the outgoing light 
ray from the outer edge of the nonzero initial data, 
i.e. from T — and rf, — a + reaches null infin- 
ity, and has to return to zero again after the ingoing 
ray from the same point, bouncing through the cen- 
ter reaches null infinity as well. This means that the 
si gnal can be nonzero only within t he time interval 
^/Kprf+l — KTb < T < ^ K^r\ + 1 + Krf, which, af- 
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FIG. 6: Logarithmic plot of the absolute error of the 
value of the massive Klein-Gordon field at constant ra- 
dius R = 0.02542 calculated with several different numer- 
ical resolutions corresponding to the indicated number of 
spatial grid points. The correct absolute value of the field 
is also plotted in order to show the time intervals where 
the error is smaller than the actual function value. We 
note that the downward pointing peaks indicate moments 
of time where the functions change signature. In order to 
reduce the complexity of the figure, for each resolution 
the error is plotted only up to a time where it becomes 
close to the exact value. 
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FIG. 7: Time dependence of the absolute value of the 
rescaled field variable z = r$ at null infinity R = 1 
obtained by our numerical code in case of the mass- 
less Klein-Gordon field with numerical resolutions corre- 
sponding to the indicated number of spatial grid points. 



ter substituting the values of the constants k and r^, 
corresponds to the interval 0.8402 < T < 1.1902. We 
see that the higher the resolution is the better the 
numerical simulation can follow the abrupt changes 
in the magnitude of the field variable. Other feature 
that can be seen from the figure is that after the sig- 
nal left the system a numerical noise remains, with 
a magnitude of about 10~^® part of the size of the 
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initial data. The size of this noise remains in this low 
range for very long time intervals, even for T « 10^. 

In order to keep the numerical noise at a low level 
for very long time intervals a smooth cutoff in the arti- 
ficial dissipation term has been introduced in the un- 
physical region R> 1. This was necessary because of 
the incompatibility of the sixth derivative dissipative 
term with the outgoing boundary condition applied 
at the outer edge of the grid. 



VII. NUMERICAL RESULTS FOR THE 
TIME EVOLUTION OF MAGNETIC 
MONOPOLES 

A. Choice of initial data 

The numerical simulations presented in this paper 
were carried out in the A = case, when the Higgs 
field is massless and the static monopole solution is 
given by Wg and Hg in fM)) and (|27l) . For our nu- 
merical simulations we used initial data which con- 
sisted of a concentrated spherically symmetric pulse 
superposed on the static monopole solution. Defining 
hs = r{Hs — Hoo), at T = we choose {w)o = Wg, 



(Wii)c 

and 



IWT 



Dr-Ws, {h)o = hs, {hR)o = OrHs, (/it)o = 



if r G [a — &, a + 6]; 



, 



otherwise. 



(91) 

In all of the following simulations we take a — 2, b = 
1.5 and d = 10, but we choose different amplitudes 
c. As it is indicated by Table. HI for instance, for our 
basic choice c = 70, which is the case discussed in (l^ . 
the energy provided by the pulse is about 55% of the 
energy of the background static monopole solution. 
Notice, however, that for c = 280 the energy of the 



c 




Ec/Es 


0.7 


0.000696595 


0.00005543 


7 


0.0696595 


0.005543 


70 


6.96595 


0.5543 


280 


111.455 


8.869 


1120 


1783.28 


141.9 



TABLE I: The value of the energy Ec of the exciting 
pulse, along with the value of its ratio Ec/Es with respect 
to the energy, Es = 12.56637, of the original static BPS 
solution is given for several values of the parameter c. 
In fact, the listed values of c are labeling exactly those 
excitations for which the properties of the time evolution 
are going to be discussed below. 

pulse gets to be 16-times of the energy of i?c=7o which 
implies that the energy of the pulse is 887% of the 
energy of the original static BPS monopole. Clearly 
such an excitation is far too strong to be considered 



as a perturbation and, as we shall see below, highly 
non-linear aspects of the evolution show up in this 
case. 



B. The qualitative picture 

Before presenting an accurate and comprehensive 
quantitative analysis of the dynamics of the excited 
monopole it is informative to consider the qualitative 
aspects of the investigated processes. This is done in 
this section by making use of spacetime plots of var- 
ious quantities, where the term 'spacetime plot' is to 
indicate that the quantity in question is given as a 
function of space and time always over a rectangu- 
lar coordinate domain of the (i?, T) plane. Accord- 
ingly, on these spacetime plots the horizontal lines 
will always represent the T — const coordinate lines, 
along which the value of the i?-coordinate is changing 
between the indicated minimal and maximal values, 
Rmin < R < Rmax, increasing from the left to the 
right. Similarly, the vertical lines are the R = const 
coordinate lines, along which the value of the T- 
coordinate is changing between the indicated initial 
and final values, Ti < T < Tf, increasing from the 
bottom to the top. 



1. The evolution of the basic variables 

This short subsection is to present some spacetime 
plots informing about the behavior of the basic vari- 
ables w and h as the monopole react for weak and 
strong excitations. 




FIG. 8: The spacetime plot of the difference w — wq between 
the dynamical and the static values of the field variable w is 
given for the weak excitation with c = 0.7 while T and -R take 
values from the intervals < T < 2.9661 and < -R < 1. The 
maximum, the minimum and the average values oiw — wo are 
0.004230, -0.002398 and 9.021 ■ 10"^. 
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FIG. 9: The spacetime plot of the difference h — ho between 
the dynamical and the static values of the field variable h is 
given for the weak excitation with c = 0.7 while T and R take 
values from the intervals < T < 2.9661 and < iJ < 1. The 
maximum, the minimum and the average values oi h — ho are 
0.0009632, -0.004673 and -0.0001738. 




FIG. 10: The spacetime plot of w is given for the strong 
excitation with c = 280 while T and R take values from the 
intervals < T < 2.9661 and < R < 1. The maximum, the 
minimum and the average values of w are 2.070, —1.210 and 
0.07346. 

There are immediate similarities and differences to 
be observed. First of all, the basic features indicated 
by Figs. [8] and [9] and Figs.flO] and [TT] are similar al- 
though the energy of the exciting pulse relative to the 
energy of the initially static BPS monopole is negligi- 
ble for the weak excitation - that is why not the field 
values w and h themselves but their difference with 
respect to the static background are plotted on Figs. [8] 
and [9] - while it is about nine times of the energy of 
the BPS monopole for the strong excitation. In spite 
of the immediate similarities the frequencies of the 
developing oscillations at the center are different. It 
is a bit of surprise that the larger the excitation is the 
lower is the initial frequency. Comparing the w and h 
evolutions, the most important differences get mani- 



FIG. 11: The spacetime plot of h is given for the strong 
excitation with c = 280 while T and R take values from the 
intervals < T < 2.9661 and < < 1. The maximum, the 
minimum and the average values of h are —5.317, ,—3.034. 

fested if one considers the shape of the curved lines 
associated with the maximum and minimum field val- 
ues. While these curves seem to approximate ingoing 
null geodesies in case of the massive Yang-Mills field, 
see Figs. [8] and \W[ they closely follow outgoing null 
rays in case of the massless Higgs field as indicated by 
Figs.[9]and[lTJ Since all physically interesting quanti- 
ties like the energy, energy current or magnetic charge 
densities are derived from w and h, along with their 
derivatives, these derived quantities shall show a mix- 
ture of these clear characters as it will be seen below. 

These figures make it also transparent that the two 
basic variables w and h, although, are strongly cou- 
pled in the central region, they can be considered as 
being independent far away from there. 



S. The spacetime dependence of energy and energy 
transfer 

One might expect that the most natural quantity to 
be applied when one tries to present the evolution of 
a system is the spacetime plot of the energy density. 
Such a quantity has a meaning only after the fam- 
ily of observers measuring it has been specified. In 
our setting it is natural to consider the energy den- 
sity measured by static observers represented by the 
unit timelike Killing vector field t"' — {d/dtY . Let us 
denote by e the associated energy density which can 
be given by the relation e = Tabf^t^ ■ 

Fig. SSI is included to make it clear that even for 
a large excitation the information concerning the dy- 
namics of the system is not properly refiected by the 
energy density. There is a clearly visible huge pulse 
at the beginning, the top of the peak is pointing far 
out from the plane of the figure at the lower left cor- 
ner, which is followed by a moderated oscillation at 
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the very center. Although the range of the radial co- 
ordinate is now only < i? < 0.4915, corresponding 
to < r < 25.92, it is easy to agree that beside the 
short dynamical starting nothing interesting is visible 
on most of this figure. 





FIG. 13: The spacetime plot oi £ for the large excitation with 
c = 280 is shown for the space and time intervals < -R < 1 
and < T < 2.9661. The maximum, the minimum and the 
average values of £ are 5572, 0, 96.68. 




FIG. 12: The spacetime plot of the energy density e = 
Tab^'^i^ is shown for a strong excitation with amplitude c = 
280 and for the space and time intervals < R < 0.4915 and 
< T < 2.119. The maximum, the minimum and the average 
values of e are 13.33, 2.214 ■ 10-^ and 0.03845. 

Investigating the system it turned out that instead 
of monitoring the behavior of the energy density e it 
is more informative to consider the time evolution of 
the energy density £ = £{T, R) 'associated with shells 
of radius R\ defined by 

dr 

f = r% = 27rAi(i?2 + l)_T% , (92) 



where the components of the energy momentum ten- 
sor ^ are calculated in the coordinate system a;° = 
(T,R,0,(p). Justification of this definition will be 
given in subsection lVII C 2} where the energy conser- 
vation will be considered in details. The integral, 

£dR with respect to the i?-coordinate gives, for 
any fixed value of T, the total energy E{T) of the 
T = const hypersurface. 

The advantages related to the use of £, instead of e, 
are visible on Fig.[l3l This spacetime plot of £ is al- 
ready very informative, it clearly manifests the most 
important characteristic features of the beginning of 
the evolution. First there is a direct energy transport 
to (represented by the right vertical edge) by the 
massless Higgs field following outgoing null geodesies. 
The second phenomenon is that right after a short 
period a quasi-stable oscillation starts at the central 
region (close to the left vertical edge). Finally, the 
formation of the expanding shells of high frequency 



oscillations in the distant region is also clearly vis- 
ible. The time dependence of the mean value, the 
amplitude and frequency of some of the oscillating 
quantities will be investigated in details in subsec- 
tion lVII C 61 

To understand the evolution of the investigated sys- 
tem it is also important to consider spacetime plots 
giving information about the energy transfers repre- 
sented by the energy current density 'associated with 
shells of radius R\ 

dr i?^ 
S = W —T\ = 27rK{R^ + i)—T\. (93) 

Referring to the results contained by subsec- 
tion lVII C2| it can also be justified that S is the 
quantity, the time integral SdT of which, for any 
fixed value of R, gives the energy passing though 
the R = const hypersurface during the time interval 
Ti<T < Ta. 

Fig.HH shows the spacetime plot of S for the case 
of the intermediate energy exciting pulse with c = 70. 
The same characteristics of the evolution as on the 
previous figure, relevant for the large exciting pulse, 
show up clearly on this figure, as well. 

In virtue of Fig.[l4]at least two different phenomena 
clearly manifest themselves. Firstly, the amplitude of 
the energy current density is much smaller beyond 
the shells of high frequency oscillations (this region is 
represented by the relatively smooth part close to the 
right vertical edge) than before reaching these shells. 
Actually, the 'ribs', representing the maximal values 
of S which apparently follow the shape of outgoing 
null rays all the way out to but they are getting 
to be modulated more and more by the slowly moving 
massive shells of oscillations in the distance. Secondly, 
it is also apparent that some sort of inward directed 
energy transfer starts at the location where the 'ribs' 



17 




FIG. 14: The spacetime plot of the quantity S, for the in- 
termediate excitation with amplitude c = 70 is shown for the 
space and time intervals < i? < 1 and < T < 2.88559. 

meet the massive shells of oscillations. 

These shells of oscillations can be nicely pictured by 
focusing to a small subsection of the previous space- 
time diagram as shown on Fig.fTSl 




FIG. 15: The spacetime plot of the quantity S, for the in- 
termediate excitation with c = 70 is shown for the space and 
time intervals 0.5829 < R < 1 and 2.11864 < T < 2.88559. 

Finally, we would like to attract attention to the 
appearance of some interesting features of the os- 
cillations close to the origin. As the strength of 
the appHed excitation is increased, the associated 
shapes, which can be clearly recognized on the fol- 
lowing spacetime plots, are getting more and more 
complex. In virtue of Figs.[l6l[l7|and[l8l it is tempt- 
ing to say that there might be a mean of switching on 
more and more oscillating degrees of freedom of the 
monopole as the energy of the exciting pulse is in- 
creased. Notice that despite the relatively late time 
interval the oscillating shapes are still present. Since 




FIG. 16: The space and time dependence of the E — So for 
the intermediate excitation with amplitude c = 70 is shown 
for the coordinate domain given by the relations < iJ < 0.35 
and 1.695 < T < 4.237. 




FIG. 17: The spacetime plot of £■ — €o, for the large excita- 
tion with amplitude c = 280 is shown for the space and time 
intervals 0<R< 0.35 and 1.695 < T < 4.237. The maximum 
values are about ten times larger than the those are on Fig. ll6l 

£ = £a ~ Q aX the left vertical edge corresponding to 
i? = 0, it can be clearly seen from the figures that for 
the larger excitations the energy contained in the cen- 
tral region is definitely smaller than the energy con- 
tained in the same region of the static monopole. This 
means that at the centre still there is a huge deficit 
in the energy, i.e. considerably large fraction of the 
energy of the static BPS monopole is still somewhere 
in the outer region. A very large pulse in this sense 
almost destroys the monopole, scattering its energy 
to faraway regions. Of course, because of the charge 
conservation this energy has to come back slowly to 
the center, forming the well locaHzed monopole again 
with smaller and smaller amplitude oscillations. 



3. The space and time dependence of the pressure 

There are two types of pressure, the radial and 
tangential, that can be associated with our spheri- 
cally symmetric system. They are given as Prad = 
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FIG. 18: The spacetime plot of £ — £o, for the extra large 
excitation with amplitude c = 1120 is shown for < -R < 0.35 
and for a later time interval 36.441 < T < 38.983. Here the 
maximum values are still more than three times larger than 
those on Fig. ll6l 



T'^'^nr'nf' and Pang = T^^nrn^;', respectively, 

where n^J^'^ and nf^ denote the unit norm space- 
like vector proportional to the coordinate differen- 
tials {dR)a and (d0)a, respectively. The radial and 
tangential pressures are not equal to each other as it 
may happen in case of generic spherically symmetric 
configurations. Note that both of these pressures are 
identically zero for the static BPS configuration. 

To be able to show the tiny features far form the 
center on these spacetime plots, on Figs.[l9] and EOl 
we needed to multiply the pressures with an enormous 
scaling factor, and use an almost orthogonal projec- 
tion (i.e. very faraway point of view) to compensate 
the scaling factor. Actually, the ratios of the values of 
the radial and angular pressure at the origin and at 
the middle of the plot, at i? = 0.5, are about 1.5 • 10^ 
and 2.3-10^, respectively. This means that both of the 
pressures are significantly larger at the origin, where 
the monopole lives, than anywhere else. 

Although there is a noticeable phase difference be- 
tween the almost sinusoidal oscillations of the two 
pressures at the origin, moreover, clear differences also 
show up in the shape of the oscillations further away, 
no really surprisingly new phenomenon, in addition to 
the formerly reported ones, can be observed on these 
figures. 



C. The quantitative picture 

Before presenting a more detailed description of the 
behavior of the fundamental physical quantities let us 
show some simple but convincing results concerning 
the numerical preservation of the constraint equations 
during the investigated long term evolutions. 




FIG. 19: The time evolution of the radial pressure, Prad, 
is shown for the intermediate exciting pulse with amplitude 
c = 70. The coordinate domain is given by the relations < 
R<1 and 16.949 < T < 19.915. Notice that this time period 
refers to that part of the evolution when the system behaves 
already much like a quasi-breather. 




FIG. 20: The time evolution of the angular pressure, Pang, 
is shown. The system and the coordinate domain is exactly 
the same as on Fig.ll9l 

1. Monitoring of the constraints 

Recall that the i?-derivatives of w and h are used 
as dependent variables, in consequence of which, the 
relations wr = daw and hn = dah play the role of 
constraints now. In the analytic setting these rela- 
tions are preserved by the evolution equations, pro- 
vided that they hold on the initial surface. One of the 
possible tests of our numerical code is the monitoring 
of the violation of these constraints. On Fig.[2T] the 
L? norm of the h constraint is presented as a function 
of time, defined as 

\\hR - dRh\\ = (^J^ {liR - dRhf dR^ ' . (94) 
It can be seen from the figure that for lower resolu- 
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FIG. 21: Time dependence of the ft constraint violation 
for several numerical resolutions corresponding to the in- 
dicated number of spatial grid points. The amplitude 
parameter in the initial data was chosen to be c = 70. 
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PIG. 22: The boundary dV = dVi U 91^2 U dVa of the shaded 
spacetime domain, represented by the union of the T = Ti , 
T = T2 and R = R hypersurfaces is shown. 



tions the error decreases when one increases the num- 
ber of grid points. However, at higher resolutions, 
most Hkely because of rounding errors, the constraint 
violation becomes higher again. It is also apparent 
that apart from a short initial period the constraint 
violation is a decreasing function of time. The time 
dependence of the w constraint violation is similar to 
that of h on Fig.[2H only its magnitude is somewhat 
smaller. 

Finally, we would like to emphasize that the numer- 
ical preservation of the constraint equations, reported 
above, in our rather long time evolutions is consider- 
ably good. Our code seem to be free, at least in case of 
the investigated system, from the deficiencies of other 
approaches, where in case of various dynamical sys- 
tems exponential increase of the numerical violation 
of the constraint equations have been observed. 

2. Energy balances 

To monitor the appropriateness of the applied nu- 
merical scheme we calculated energy balances for var- 
ious spacetime domains. In particular, we considered 
the type of domains shown by Fig. [221 bounded by 
T = const, and R = const, hypersurfaces. 

To see how the energy balances can be calculated 
recall that whenever t"" is a Killing vector field, satis- 
fying the Killing equation V^'^f^-' = 0, the vector field 

f = T\t'' (95) 

is divergence free, i.e. 

Vaf = . (96) 
Then, by making use of Stokes' theorem we have 

/ Vafe= f na.fi = 0, (97) 

Jv JdV 



where e denotes the 4-volume element while e is the 
3-volume element induced on the boundary dV of V. 
Note that iabc — nee'^abc where n^ is the (outward 
pointing) unit normal 1-form field on dV. For the fol- 
lowing considerations we choose t° as the unit normed 
static timelike Killing vector field t°- = [d/dty. 

Let us denote by B{T, R) the ball of radius R = R 
centered at the origin of the T = T{= const) hyper- 
surface, moreover, by C(T2,Ti,_R) the portion of the 
cylindrical hypersurface R = R between the T — Ti 
and T — T2 hypersurfaces. With these notations the 
boundaries, as indicated on Fig. [221 dVi and dVs are 
the balls B{Ti,R) and B{T2, R), while dV2 is the cylin- 
der C{T2,Ti,R) connecting them. Then the energy 
contained in a ball S(T, R) can be given as 

E{T,R)= f n^:^fi^'\ (98) 

J8{T,R) 

where e^*-' is the volume element on the constant T 
hypersurface and ni* is its future pointing normal 
vector. Similarly, the energy transported through the 
timelike boundary C{T2,Ti, R) is given as 

S{T2,Ti,R)^ [ ni'^r~e^'\ (99) 

Jc{T2,Ti,R) 

where e^*^ is the volume element on the constant R 

(s) 

hypersurface and na is its outwards pointing nor- 
mal vector. Using these notations the energy balance 
equation (|97l) takes the form 

E{T2,R)-E{Ti,R) + S{T2,Ti,R)=0 . (100) 

Note that while in analytic considerations the l.h.s. 
of (|100p is always identically zero, for any choice of 
T2, Ti and R, in meaningful numerical simulations the 
corresponding quantity is merely close to zero, i.e. 
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an apparent violation of the energy conservation hap- 
pens. In fact, the numerical value NyeciT2,Ti, R) of 
the "violation of the energy conservation" associated 
with a spacetime domain V, which is defined as 

N^ec = E{T2, R) - E{Ti,R) + S{T2, Ti,R), (101) 

can be used as one of the possible monitorings of the 
appropriateness of a numerical code. 

To evaluate the integrals l|98p and ijOQl) we need to 
determine the volume elements e*^*^ and e'^-* which can 
be given as special cases of the relation 



^abc 



g{4)\£eabcn 



(102) 



where Eeobc denotes the Levi-Civita alternating tensor 
with £0123 = 1- To get e^*^ and e*^^^ consider now 
the Minkowski spacetime in coordinates (T, i?, 0, 4>) 
associated with the appHed conformal representation. 
Then the components gap of the metric tensor, gab, 
and its inverse read as 
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(103) 



and 
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" fl2 sin^ e / 

(104) 

Calculating the determinant 3(4) of the spacetime 
metric we obtain 



|5(4)| 



l)i?2 sin 61. 



(105) 



The future pointing normal form of S(T, R) has 
the components n^a = (l/\/ g'^^, 0, 0, 0). This implies 
then that n(*) ° = g°°n~Q^ = 
we have that 
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9 = W+i-^ moreover. 



-it) 



i?2 sin 6* 



(106) 

where h^-*' is the determinant of the three metric 
i-j induced by gab on B(T, R), i.e. 

i?2 sin 6* 



'*ab — 9ab — na 



\hW\ 



(107) 



on B{T,R). Similarly, the normal of C(T2,Ti,i?), 

(s) 

pointing out from the domain V, is no — 
(0,l/v^-5ii,0,0) and thereby n(^) ^ = g^'^n'f^ = 
^\f--g^ = — which implies, in particular, 

that 



i?2 sin 6* 



J|/i(^)|(dr)„A(d0)^A(d</))^, 

(108) 



where ft-^"-' is the determinant of the induced metric 

i?2 sin 6* 



"afc ^ 5ah + "•a «b , and 



(109) 



Taking the above relations into account we get that 

E{T,R)^ £{T,R')dR' , (110) 
Jo 

S{T2,Ti,R)= [ 's{T,R)dT , (111) 



and 



where £ and 5 are the functions defined in Ip2|) 

On Figs. [23] and [24] the time dependence of the 
violation of the energy conservation Nyec{0,T, R) is 
shown for numerical runs with various resolutions 
evolving from initial data with amplitude c = 70. 
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FIG. 23: Initial behavior of the numerical violation of 
the energy conservation Nvec(T, 0, R) for different spatial 
resolutions is shown. The spacetime domain is bounded 
by the constant time surfaces Ti = and T2 = T, while 
the cylindrical boundary is at 7? = 0.83051, correspond- 
ing to r = 107.07. 

It can be seen from these figures that the error in the 
energy conservation decreases according to expected 
fourth order convergence during an initial period up 
to approximately T = 0.7 and also later, after T > 20. 
In the intermediate region the conservation violation 
decreases more slowly. This is possibly due to the 
highly oscillating character of the fields in that time 
interval, especially since the energy current integrals 
were calculated only with a second order convergent 
method. 

The magnitude of the presented constraint viola- 
tions is really meaningful only if it can be compared 
to the full energy content. On Fig.[25l the time de- 
pendence of the energy content E{T, R) and energy 
loss S{T, 0, R) is shown for the same R as on Figs. [23] 
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FIG. 24: The time dependence of the violation of the 
energy conservation Nvec{T, 0,R = 0.83051) for the same 
numerical simulations as on the previous figure is shown 
but for a much longer time interval. 

and [241 To bring the two curves of the graph into the 
the same range, instead of the total energy E{T, R) 
the energy difference between E{T, R) and that of the 
static monopole solution Es{R) is shown on Figs. [25] 
and[26l For the chosen radius R — 0.83051 the energy 
of the static background monopole is Es{R) — 12.449, 
less than double of the energy provided by the initial 
deformation. 
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FIG. 25: The solid curve on the graph shows the time 
dependence of the energy contained inside a radius R on 
the constant T hypersurface E{T, R) decreased by the en- 
ergy content of the static monopole, Es{R). The dashed 
curve represents the total energy transported through 
the constant R hypersurface until time T. The sum of 
these two curves is constant up to errors represented by 

iV„ec(T,0,i?). 

The time dependence of the difference E{T, R) — 
Es{R) of the dynamical and static energy functions, 
E{T,R) and Es{R), for longer time periods is shown 
logarithmically on Fig.[26| for increasing numerical 



resolutions. The convergence of the presented curves 
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FIG. 26: Time dependence of the extra energy content 
E{T, R) — Es (R) for the same system as on the previous 
figure, but now for longer time periods and with logarith- 
mic axes. 

provides us a strong indication on how long we can 
take our numerical results seriously. For the highest 
resolution used, i.e. that with 2^^ spatial grid points, 
the calculation becomes unreliable after T — 300. We 
note that in physical time, i.e. in time measured in 
mass unites, the time interval AT = 300, at the cen- 
tre, corresponds to At = 6000, meanwhile the central 
monopole performs nearly one thousand oscillations. 



3. Linear and non-linear effects 

In order to distinguish linear and nonlinear ef- 
fects in the evolution of magnetic monopoles we first 
present results corresponding to very small initial 
deformation, with amplitude c — 0.7. Then the 
extra energy provided by the initial pulse is only 
0.000696595, which is very small compared to the 
energy of the static monopole, which is 12.56637 to 
seven digits precision. On Fig. [27] the time depen- 
dence of the extra energy is presented in a ball of ra- 
dius R — 0.067797, which corresponds to r = 2.7244. 
Inside this radius the static background solution con- 
tains energy Es{R) — 7.0238, which is about 56 per- 
cent of the total energy of the static monopole. For 
this low energy case, after a short initial period, the 
energy content starts oscillating around the energy 
of the static monopole. The time decay of these 
oscillations is r~5/6 ^ good approximation. On 
Fig. [28] the energy content is presented in a much 
larger ball, with radius R — 0.83051, corresponding 
to r = 107.07. Inside this large radius the energy con- 
tained in the static monopole is 12.449, which is 99 
percent of the total static energy. It is apparent that, 
because the Yang-Mills field is massive, part of the en- 
ergy provided by the pulse remains inside this radius 
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r =2.7244 
Es(R) =7.0238 
c = 0.7 , E(0,R) - Es(R) = 0.00069208 
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E(0,R)-Es(R) = 0.00059559 




FIG. 27: The extra energy E{T, R)-E,{R) as a function 
of time contained inside radius R = 0.067797 is shown 
for an initial and for a longer time period, 1 < T < 300. 
By a careful analysis it can be justified that on average 
the energy content close to the origin might be slightly 
lower than that of the static monopole. The forming of 
a long lasting breathing phase with a clear asymptotic 
time dependence can also be clearly recognized. 



FIG. 28: The extra energy content E{T,R) - Es{R) 
for radius R = 0.83051 is presented in case of the low 
energy deformation, with c — 0.7. It is transparent that, 
although the energy content close to the origin might be 
slightly lower than that of the static monopole, up to this 
relatively large radius the system possesses more energy 
on the average than the static monopole had. 



for quite a long time. The oscillations of the energy 
content will be centered around the static value only 
after about T > 100. 

On Fig. [29] the energy contained in the same small 
radius R = 0.067797 as the one used at Fig. [27] is 
shown for a large initial deformation corresponding 
to c = 280. In this case the energy provided by the 
pulse inside this radius is 110.73, which is about 15 
times more than the energy of the static monopole 
in this region. The main difference compared to the 
previous case is that the oscillation of the energy is 
not centered on the static monopole energy anymore. 
The average value of the energy is below the static 
value for a very long time period. We interpret this as 
a nonlinear effect. The expanding initial pulse sweeps 
out not only the provided extra energy but also a part 
of the energy of the static monopole. It takes a long 
time to get back this energy which is stored in the 
massive Yang-Mills field oscillations at intermediate 
distances from the monopole. On Fig. [30] the extra 
energy is shown for the same large initial pulse but 
for the larger radius, r = 107.07. Inside this large 
radius, just like in case of the tiny excitation with c = 
0.7 (see Fig. [27j) . the energy remains above the static 
value during the entire evolution where our numerical 
simulation can be considered to be vaHd. 

On Fig. [31] the effect of an extra large initial pulse 



is presented. In this case the amplitude of the initial 
pulse was chosen to be c = 1120, corresponding to 
provided energy 1771.75 inside the R < 0.067797 re- 
gion. In this case the initial pulse quickly sweeps out 
most of the energy of the monopole from the central 
region. However, after some highly energetic oscilla- 
tions, the major part of the monopole energy returns, 
and the energy content starts to oscillate around a 
value below that of the static monopole, like in the 
case of the less energetic nonhnear excitations consid- 
ered previously. 



4- The time dependence of the total radiated energy at 

In this subsection we provide a short account on the 
time dependence of the energy radiated to future null 
infinity. Since the Yang-Mills field is massive, w de- 
cays exponentially when approaching the i? = 1 line 
corresponding to null infinity. This implies that the 
intensity of the radiation to J^"*" is determined only 
by h and its first derivatives at i? = 1. The full radi- 
ated energy up to the moment T, i.e. St = S{T, 0, 1) 
in terms of the notation introduced in section lVII C 2\ 
depends on the amount of energy Ec provided by the 
initial deformation ijOTI) . For our special choice of ini- 
tial data Ec is proportional to c^. On Fig.[32]the ratio 
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FIG. 29: The extra energy content E{T,R) - Es{R) 
for radius R = 0.067797 is shown for high energy de- 
formation corresponding to c = 280. The graphs clearly 
indicate that for a large exciting pulse the average energy 
content of the central region is below that of the static 
monopole. 



FIG. 31: The full energy content E{T,R) inside radius 
R = 0.067797 is indicated for an enormously large initial 
pulse with c = 1120. Even though the monopole seems 
to be swept out from the central region for the initial 
period eventually the behavior of the system returns to 
an oscillating normal monopole state. 
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FIG. 30: The extra energy content E{T,R) - Es{R) is 
shown for the larger radius, R — 0.83051, and for the 
high energy deformation with c = 280. The asymptotic 
power law time decay can easily be read off the log scale 
plot. 

St/Ec is plotted for five different initial data ampli- 
tudes. It can be seen that except for very large initial 
energies the radiated energy is closely proportional to 
the energy contained in the initial pulse. In Table ITTl 
for the five chosen c amplitude value, we list the en- 
ergy of the initial pulse Ec, its ratio to the energy of 
the static monopole, Ec/Eg, the fraction of the ini- 
tial energy radiated up to T = 300, and the fraction 



of the initial energy radiated during the whole time 
evolution. The value ^oo was estimated by assuming 



c 


Ec 


EcjEs 


Ssoo/Ec 


Soo/Ec 


0.7 


0.000696595 


0.00005543 


0.47793 


0.48014 


7 


0.0696595 


0.005543 


0.47428 


0.47646 


70 


6.96595 


0.5543 


0.44037 


0.44250 


280 


111.455 


8.869 


0.56388 


0.5725 


1120 


1783.28 


141.9 


0.85177 


0.8555 



TABLE II: The energy Ec of the initial pulse, its ratio to 
the static monopole energy, i.e. Ec/Es, furthermore, the 
fraction of the energy radiated to J^'^ up to time T — 300 
and T — oo, are given for five different values of the initial 
amplitude parameter c. 

that the energy content of the system decreases ac- 
cording to a power law decay with T~^^^ during the 
later stages of the radiation process. As it can be 
seen on Figs. [33] and [Ml this is certainly a very good 
approximation for weaker initial excitations, and it 
appears to be consistent with the late time behavior 
of more energetic evolutions. 

Notice that since < Ec the limit of the energy 
associated with the system, on T = const hypersur- 
faces, does not tend to the energy of the original BPS 
monopole. This observation is in fact in accordance 
with the existence of the massive shells of high fre- 
quency oscillations. These shells never reach future 
null infinity, and hence they store for the rest of the 



24 



O.f 
0.7 
0.6 
0.5 
0.4 
0.3 













0=0.7 




0=7 

0=70 




0=280 
c=1120 



0.5 



1.5 

T 



2.5 



0.9 
0.8 
0.7 
0.6 
0.5 
0.4 
0.3 
0.2 
0.1 




0=0.7 
0=7 
0=70 
0=280 
0=1120 



50 



100 



150 
T 



200 



250 



300 



FIG. 32: The ratio St/Ec is plotted for five different 
choices of the initial amplitude c, where St ~ S{T, 0, 1) 
is the energy radiated to infinity up to T and Ec denotes 
the energy provided by the exciting pulse. 
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FIG. 34: The quantity (Soo -5't)/£c multiplied by T^/^ 
is plotted for the five different initial data configurations. 



5. The space and time dependence of the energy 
transfer towards J''^ 
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FIG. 33: The quantity {So^-St)IE^ is plotted logarith- 
mically for the indicated values of the initial amplitude 
c. 

entire evolution a part of the energy of the original ex- 
citing pulse. This phenomenon gets to be even more 
transparent if one thinks of a simple massive Klein- 
Gordon system starting, say, from a trivial configura- 
tion where only the exciting pulse stores energy. The 
pulse disperses, as it is expected, but the system is 
conservative so the total energy has to be conserved. 
Due to this, there is no upper bound on the increase 
of the frequency associated with these shells of oscil- 
lations forming even in this simplest possible linear 
case. Let us mention that a more detailed discussion 
of this phenomenon can be found in . 



One of the possible methods in monitoring the 
main features of the energy transfers is to investi- 
gate in more detail the space and time dependence of 
the extra energy content of the dynamical monopole. 
Fig. [35] pictures the time dependence of the difference 
A£: = E{T,R) - Es{R) along the constant i?-slice, 
with R ~ 0.1949, corresponding to r = 8.1045. In 
particular, on the upper part of the figure the time 
dependence of AE, along with its maximum and min- 
imum contours, together with the associated aver- 
age or mean value, are plotted for the time inter- 
val 2 < r < 17. Notice that, even in case of the 
considered intermediate excitation, at the indicated 
relatively long time period the mean value of the dif- 
ference E{T, R) — Es{R) is still negative. Looking at 
this plot it is also clear that the time dependence of 
the amplitude of the high frequency oscillations and 
that of the mean value should be investigated sepa- 
rately. For this purpose we investigated separately 
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FIG. 35: On the upper plot the difference AS = 
E{T,R) — Ea{R), along with the maximum, minimum 
contours and the mean values, is plotted for the time in- 
terval 2 < T < 17. On the lower plot the maximum, 
minimum contours associated with the separated high 
frequency oscillations, i.e., with AE — {AE), are shown, 
where {AE) denotes the mean value of the difference 
AE = E(T, R) - Es{R). 



the time dependence of the mean value {AE), and 
the time dependence of the amplitude of the high fre- 
quency oscillations (AE' — {AE))max- On the lower 
part of Fig. [35] the behavior of these separated vari- 
ables are shown for the same time period as on the 
upper graph. 

By making use of log-log plots, see e.g. Fig.[36l for 
various i?-slices it is again straightforward to recog- 
nize a power law time decay that can be associated 
with both the amplitude of the high frequency oscilla- 
tions, {AE — {AE))max, and that of the mean value, 
{AE). In particular, denoting by —a and the two 
power law exponents we have that {AE) ^ and 
{AE - {AE))^ax ^ T-", respectively. 

Table [ini collects more information about the in- 
dependence of the values of the power law decaying 
exponents a and in case of excitation with c = 70. 
Apparently there is a certain universality in the dy- 
ing out of the high frequency oscillations, where the 
numerical value of a is very close to 5/6. Neverthe- 
less, contrary to our expectations, it turned out that 
the time dependence of the mean value {AE) reflects 
some sort of non-universal features, i.e., the value of ^ 
varies, so it is, in general, different from 1.5 which had 
been associated with the R = 0.8305-slices in case of 
the excitations with amplitudes c = 70 and c = 280 
(see Figs.EiandlSO]). 



r 


R 


M 


a 


1.0176 


0.0254 


0.9 


0.83333 


2.7244 


0.0678 


0.94 


0.83333 


4.4609 


0.110169 


0.9 


0.828 


6.2471 


0.1525 


0.87 


0.833 


8.1045 


0.1949 


0.86 


0.83 


13.232 


0.3010 


1.2 


0.8333 


19.497 


0.4068 


1.5 


0.83333 


27.822 


0.5127 


1.52 


0.833 


40.088 


0.6186 


1.54 


0.8333 


61.018 


0.7246 


1.55 


0.833 


107.074 


0.8305 


1.53 


0.84 


146.655 


0.8729 


1.5 


0.877 



TABLE III: The asymptotic time dependence of {AE) 
and {AE — {AE))max, along R = const slices, were found 
to follow the power law time decay {AE) ~ T~'^ and 
{AE — {AE))max ~ T~", respectively. This table provides 
the values of n and a for the indicated 7?-slices in case of 
excitation with c = 70. 



It is important to emphasize that by the mere na- 
ture of the applied approximations there is no nu- 
merical method which could represent the full history 
of the expanding shells of high frequency oscillations 
properly. Therefore it is of distinguished importance 
to know how strongly the above reported results, con- 
cerning the energy transport, depend on the appHed 
numerical resolution. Is there at all a meaningful con- 
vergence guaranteed, say, at a suitably large value of 
R close to the region occupied by the shells of high 
frequency oscillations? Fig.[36]is to justify that for the 
largest value of R referred to in Table lllH i.e. at the 
R — 0.8729 slice, such a convergence is guaranteed by 
our numerical method. 
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FIG. 36: The convergence in the log-log plots of AE 
and {AE — {AE))max is demonstrated along the i?-slice 
with R = 0.8729 and for the time interval 40 < T < 500 
shown. This log-log plot clearly manifests also the power 
law decay of the amplitude of the high frequency oscil- 
lations, {AE — {AE))max, and that of the mean value, 
{AE), with the particular values a = 0.877 and /i = 1.5. 
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It is important to interpret the above findings in 
short physical terms. For instance, Table IIIII clearly 
justifies that there is an inner region with certain type 
of universal time decay for the mean value of the ex- 
tra energy content which can be characterized by the 
value fj, 0.86 — 0.9. Essentially, this region corre- 
sponds to the core of the monopole. Since the value of 
^ is the smallest in this inner region the energy trans- 
fer remains at a much higher level here for the entire 
evolution. It is also interesting that the value of ^ is 
smaller a bit apart from the centre than close to the 
centre which indicates that there remains some energy 
bouncing back and forth during the entire evolution 
in the central region where the monopole lives. Then, 
there is a wide region at the middle, characterized by 
the interval 0.3 < R < 0.87, where the mean value 
of the extra energy content is decaying more rapidly 
than elsewhere. Notice that the increase of the ex- 
ponent from about R — 0.3 till R = 0.72 with 
the increase of the value of R indicates that larger 
and larger fraction of the extra energy provided by 
the exciting pulse is stored by the expanding shells of 
high frequency oscillations of the Yang-Mills field. In 
fact, this conclusion is also supported by the slower 
and slower time decay as the value of R is further 
increased above R — 0.72. 

Let us finally investigate the energy content stored 
in expanding shells of high frequency oscillation 
formed by the massive Yang-Mills field represented 
by w at large radiuses. The energy content outside 
of R is E{T, 1) - E{T, R). Since initially the energy 
contained in this region is Es{l) — Es{R), the energy 
stored in the shells is 

-Bshoiis - E{T, l)-EiT,R)-[Esil) - ^.(i?)] (112) 

which can be written as 

^^shciis = [E{T, 1) - £;(oo, 1)] - [E{T, R) ~ E,{R)] 
+ [E{(x,,l)-Es{l)] . (113) 

Recall (see subsection lVII C 4|) that the mean value of 
the total energy E{T, 1) associated with T = const 
sHces is approaching its non-zero Hmit value E{oo, 1) 
by following a time decay of the form T~^/^. As 
it can be seen from Table IIIH the second term on 
the right hand side of (|113p always dies out faster 
than the first one. This relation justifies that the 
energy stored by the shells of high frequency oscil- 
lations has to approach asymptotically the positive 
value E{co, 1) — £'^(1) from above with the power law 
decay ~ T-^/s. 

6. The time dependence of the frequency of oscillations 

The time dependence of the frequency of the mas- 
sive w and massless h fields display markedly different 
characteristics. Frequency dependence of composite 



quantities such as the energy density will show a mix- 
ture of these properties. As can be seen on Figures 
[9] and [TH the behavior of the Higgs field h is essen- 
tially the same at each radius R, apart from a slowly 
increasing phase shift going outwards. But even this 
time delay is not very significant, since this change is 
along outgoing null geodesies. A null geodesic ema- 
nating from the center R = reach null infinity R = 1 
in a coordinate time interval AT = 1, which is small 
compared to the length of the simulations performed 
by our numerical code. The behavior of the w field 
is qualitatively different, as it is apparent from Fig- 
ures [8] and [TOl In the central region the two fields are 
strongly coupled and consequently their frequency is 
the same. On Fig. [37|the time dependence of the fre- 
quency of the w field is shown for four different initial 
data at a fixed radius R — 0.02542, corresponding to 
r = 1.017. The time dependence of the h field at this 
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FIG. 37: Time dependence of the oscillation frequency 
u} of the field w at R = 0.02542, i.e. at r = 1.017, for 
four different initial data labelled by the strength param- 
eter c. Since the frequencies approach 1 from belo-w, the 
value of 1 — a; is plotted logarithmically. The evolution 
of the frequency for the weak (c = 0.7) and intermediate 
energy (c = 70) excitations tends as r"^/^ to the limit 
value 1. For the higher energy simulations 1 — to falls 
off more quickly. For the c = 1120 simulation the down- 
■wards pointing peaks correspond to moments -where the 
frequency approaches 1 closely from below and then gets 
smaller again. 

radius would yield an essentially identical figure. 

The difference between the frequency evolutions of 
the w and h fields gets manifested at higher radii, far 
from the core of the magnetic monopole. At these 
distances expanding shell structures appear in w, os- 
cillating with high frequencies. We know in virtue of 
the results of that at a given radius R the fre- 
quency of these shells is the highest just after direct 
outgoing geodesies emanating from the initial pertur- 
bation region of the T = initial hypersurface reach 
out to the coordinate radius R. Of course, because 
of the high coordinate velocity of the outgoing null 
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rays, this happens before T = 1. It is an interesting 
result of that if R approaches 1 then the high- 
est frequency increases without any bound. After the 
highest frequency expanding shells have left the ra- 
dius R the frequency is decreasing steadily, until it 
falls to a minimum value below 1, and then starts to 
approach 1 from below, more or less according to the 
T~^/^ law observed in the central region. The time 
evolution of the frequency of the w field for an initial 
excitation with c = 70 is shown at five different radii 
on Figures [381 and [391 
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FIG. 38: Time dependence of the frequency of the os- 
cillation of the w field at five different constant radii for 
a shorter time interval, 1 < T < 7. The parameter in the 
initial data was chosen to be c = 70. Outside the core of 
the monopole the frequency decreases from a peak value 
reached before T = 1. Although the T < 1 region is 
not shown because of the high error of our frequency de- 
termination method there, the frequency reaches at least 
ij = 10 for the early stages at larger radii. 



7. The magnetic charge 

It is known that in a gauge theory, like the one 
considered in this paper, by making use of the for- 
mal definition applied in Maxwell theory, a conserved 
magnetic charge can be defined although the yielded 
charge is a Lie-algebra-valued quantity. Hence a 
meaningful definition of magnetic charge requires the 
use of a gauge independent specification of what is 
considered to be the "electromagnetic content" of the 
system. 

Historically two main approaches developed. The 
first one was suggested by 't Hooft [2^ . It is based 
on the use of the gauge independent expression [27] 

J'ab - j^tri^jFah) + -^tri^lDai;, Di,^]), (114) 

where — tr{^pip)i . It can be checked that in terms 
of the variable Aa — j^tr{il;Aa) the above "Maxwell" 
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FIG. 39: Time dependence of the frequency is shown 
for the same system and same choice of radii as on the 
previous figure, but for a longer time interval, 1 < T < 
300. Since the frequencies tend to 1, the absolute value of 
1 — ai is plotted logarithmically. The downward pointing 
peaks correspond to moments where lu falls below 1. This 
happens later for larger radii. 

tensor can be given as 

J^ab ^ daAb - dbAa- (115) 

Note, however, that in spite of the elegance of this 
construction the main drawback is that according to 
this definition the magnetic charge can reside only at 
the zeros of the Higgs field, which means that in our 
case the magnetic charge has to be point-like, it is 
concentrated to the origin. 

According to a general expectation the 't Hooft- 
Polyakov magnetic monopole has the preferable prop- 
erty of being non-singular in contrast to the Dirac 
monopole (see e.g. jlal). A proposal, made by 
Polyakov and Fadeev, fitting to this desire, which will 
also be appHed in the rest of this section, is given by 
the simple relation {3] 

^ab = -^tri^Fab) , (116) 

-no 

where Hq is the vacuum expectation value of the 
Higgs field (see Eq. 

This section is to provide a brief account on the 
associated definition of the time and location depen- 
dent magnetic charge density and the conserved total 
magnetic charge. Note that in case of the applied min- 
imal dynamical generalization of the static t'Hooft- 
Polyakov magnetic monopole no electric charge or 
electric charge density are generated. However, in 
more generic considerations, involving all the physi- 
cal degrees of freedom, associated with e.g. the most 
general spherically symmetric excitations of the BPS 
magnetic monopole, non-trivial electric charge den- 
sity appears, although the vanishing of the total elec- 
tric charge is guaranteed 
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Analogous to the arguments applied in Maxwell 
theory the magnetic charge density can be calculated 
by making use of the magnetic field strength 



Ba — — -€abcd u'' J^'^ 



(117) 



which is defined with respect to a family of observer 
represented by a timelike unit norm 4-velocity field 
m". The magnetic charge density is given then as the 
3-divergence of Ba 



„M 



(118) 



where Vq denotes the covariant derivative operator 
associated with the 3-metric induced on the 3-space 
"orthogonal" to the 4-velocity field u°-. 

To be more specific and precise, let us choose the 
timelike unit 4-vector field Ua to be everywhere nor- 
mal to the T — const, hypersurfaces. Then, for the 
components of Ua 



1 



0, 0, 



r 1 ^5 ^5 



holds, while the induced metric on the T 
hypersurfaces can be given as 

hab = gab - UaUb. 



(119) 



const. 



(120) 



Notice that for any fixed value of T the vector field 
coincides with the timelike unit norm n'-*-'" introduced 
in the previous section. 

Carrying out all the necessary calculations and 
combining all the above relations we get that 

i^.^-^^d-^^), (121) 

moreover, the expression of the magnetic charge den- 
sity simplifies to the relation 



p^^{T, R) 



gHoR^ 



dR{H[l-w^]), (122) 



where, in virtue of (|28|) and l|36p . the function H 
stands for the expression H = hQ/R+Ho. It is impor- 
tant to note that in the smooth setting the function 
dfi (^H[l — w^]) vanishes up to second order in a small 
neighborhood of the origin, i.e. there exists a smooth 
function <f = (p{T, R) so that dn^ {H[l - w^]) = LpR'^. 
This, however, guarantees that the magnetic charge 
density is regular at the center as it was antici- 
pated. 

Let us define, again referring to analogies from 
Maxwell theory, the magnetic charge associated with 
a region V in a spacelike hypersurface S by the inte- 
gral 



(V) 



(123) 



where now e denotes the 3-volume element of the 
spacelike hypersurface S. 

Returning to the previously used specific choice let 
us choose V to be the ball of radius R centered at the 
origin on a T = const, hyperboloidal hypersurface, 
i.e. V = B{T,R). Notice, that according to ifnU]) 
the components of the 3-volume element e takes the 
form given by the relation l|106p . Moreover, in virtue 
of the Stokes' theorem and (|118p the magnetic charge 
associated with the region B{T, R) can be given as a 
surface integral on the boundary of B{T, R) 

Q*^(r,i?)=/ V'^B.i^ f ly'B.e, 

Jl3(T,B.) JdB(T,B.) 

(124) 

where ly"' denotes the outward pointing spacelike unit 
normal field on the boundary dB(T, R) of the region 
B{T, R), while e is the volume element associated with 
the induced metric Xab on dB{T, R) which is given as 



Xab = hab + VaVb. 



(125) 



Since the outward pointing unit normal field v"" in 
the considered special case is proportional to (^) 
we have for the components of v"' that 



1 



= ,0,0 , 



and hence, also that 



1 f -R^ 







R^ sin^ e 



(126) 



(127) 



which leads to the relation 
i?2 sin 



ti sm a f : 

= ^2 = y'|X(2)I A [d<i))p, (128) 

where x^2) the determinant of the 2-metric Xab^ i-e. 

i?2 sin 6* 



(129) 



Combining all the above relations, for the magnetic 
charge Q*^(T, i?) associated with a ball of radius R 
and centered at the origin on a constant T hypersur- 
face we get 



4tt H , 

= il-w 

9 Ho ^ 



gHo 
2^ 



sin 9 dO < 



(130) 



The total magnetic charge, Q*^, associated with the 
entire of a constant T hyperboloidal hypersurface, can 
now be given as 



= Q*^(T,i?=l), 



(131) 
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which, in virtue of (|130p and since w tends to zero 
exponentially, while H/Hq tends to one as i? — > 1, 
takes the value 

= y . (132) 

Accordingly, the total magnetic charge Q*^ is a time 
independent quantity so that its value coincides with 
the value of the topological charge of the magnetic 
monopole. Notice also that in virtue of l|122p by excit- 
ing the the initially static BPS monopole by making 
use only of a non-zero time derivative of the Yang- 
Mills variable no magnetic charge density, and hence, 
in accordance with the above conclusion, no magnetic 
charge is added to the system. 

On Fig. HO] the spacetime plot of p^^ , the magnetic 
charge density associated with the shells of radius R, 
relevant for the intermediate initial pulse with ampli- 
tude c = 70, is shown. This quantity, /5*^, is defined 
analogously to those used previously. Accordingly, on 
any T = const hypersurface the integral /p p^^ dR is 
supposed to give the magnetic charge contained by 
the ball of radius R, and as it can be checked easily, 
-M _ ^ p'^ y/\h^\d9 dtp or, even more directly, 

= 47rr!3/-R' • = 4^/(5i?o) • dn {H [l - w^]). 
In accordance with the facts that w and H both tend 
to their asymptotic values, i.e. to zero and one, re- 
spectively, very rapidly, even in case of dynamical sit- 
uations a significant part of the total magnetic charge 
remains in a small size central region. 




FIG. 40: The space and time dependence of p*^, the mag- 
netic charge density associated with the shells of radius R, is 
plotted, for initial data with amplitude c = 70, for the space 
and time intervals 0<i?<landO<T< 4.237288135. It is 
transparent that, even though considerably large energy trans- 
ports occur during the dynamical part, the magnetic charge 
remains concentrated to the central region throughout the en- 
tire evolution. Moreover, an intensive inward pointing current 
can also be observed in the distant region which quickly re- 
turns the lost magnetic charge to the center. Accordingly, the 
amplitude of the oscillations is decreasing very rapidly there. 



This conclusion is also supported by the lower 
graph of Fig. mi where the time dependence of 
(T, R) , the magnetic charge contained by the ball 
of radius R = 0.83051 (r = 107.07) centered at the 
origin, is shown. Initially the value is constant un- 
til the arrival of the exciting pulse which sweeps out 
a small fragment of the magnetic charge. Neverthe- 
less, the missing charge returns very rapidly according 
to the two phase power law scheme indicated by the 
graph at the bottom. 
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FIG. 41: On the upper graph the time dependence of 
the magnetic charge density p^' is shown at the constant 
7?-line with R = 0.067797 (r = 2.7244). On lower graph 
the time dependence of the magnetic charge contained by 
the ball of radius 7? = 0.83051 (r = 107.07), centered at 
the origin, is shown. Here Q'^ q denotes the value of the 
magnetic charge contained by the same ball before the 
arrival of the exciting pulse. Notice that Q*^o is, in fact, 
the value of the magnetic charge contained the same ball 
in case of the static monopole. 
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Appendix A 



Consider now the expansions 



The regularity properties l|17p - (fTOl) of the field 
variables at the origin can be justified by substituting 
the expansions 

n 



fc=0 



n 

^(*'^) = E M^fcW^' + ^^(^"^')' (A-2) 



fc=0 



into the field equations and requiring that the coef- 
ficients of the various powers of r vanish identically. 
The yielded restrictions are 



wo{t) - 1 

d';wl^, ^ wk{t) = 0, 



for k being odd, and 



Hait) = 
d^Hi^„ = Hk{t) = 0, 



(A.3) 
(A.4) 

(A.5) 
(A.6) 



for k being even. It also follows from these relations, 
along with IjA.ip and (|A.2p , that at the origin the time 
derivatives dtw and dtH have to vanish throughout 
the evolution. 

Note that part of these restrictions can also be de- 
duced in a different way. In particular, by substituting 
(|A.ip and l|A.2p into the tetrad components of the en- 
ergy momentum tensor we immediately get that these 
components cannot be finite at the origin unless the 
relations 



wo(t) 

Wl{t) 

Ho{t) 



(A.7) 
(A.8) 
(A.9) 



are satisfied. 

To get the regularity result at future null infinity 
we need to transform the field equations (fTSl and (fT6| 
into a form suitable to study the asymptotic behavior 
of the fields at J^^ . Therefore, instead of the standard 
coordinates [t, r, 6*, </>) we apply coordinates (m, x, 6, cf)) 
based on outgoing null geodesic congruences deter- 
mined by the relations 



u 



t 



and X ~ —. 

r 



(A.IO) 



In this frame future null infinity is represented by the 
X = hypersurface, moreover, the field equations (fT5| 
and iflGl) read as 



x'^d^w + 2xdxW + 2dxduW — 



{w^ - 1) + 



x-'diH - 2duH -f 2xdudxH 



xH 



2w^ 



2x 



(A.ll) 



,(A.12) 



" 1 

Y,-^^k{u)x'' + 0^,{x^+^), (A.13) 



fe=0 



n 

H = Y,-Hk{u)x^ + OH{x''+^), (A.14) 

fe=0 

which are valid in a neighborhood of x = whenever 
w and H are at least C" functions through J^^ , where 
the notations Wk{u) = dxw[u,x)\x=[) and Hk{u) = 
x)|a;=o have been used. By substituting these 
expansions into the field equations l|A.lip and l|A.12p 
we get the following: Whenever X^O 



w{u) = 
H{u) = 

dtH = 





Ho 

Wk{u) = 
Hk{u) = 



(A.15) 
(A.16) 
(A.17) 
(A.18) 



for all < A; < n, i.e., in accordance with the re- 
sults of Winicour [3^1, both of the fields decay faster 
than x" to their limit values at J^"*". On the other 
hand, whenever A = and H^ is non-zero only the 
vanishing of the x-derivatives of w (up to the order 
of n) is guaranteed while we get duH = 0, i.e. H^ 
is independent of u, there is no restriction for dudxH 
and, finally, du{dxH) — for 2 < k < n provided 
that singular behavior of the field at future (or past) 
timelike infinite is excluded. Note that the missing of 
a restriction on dudxH is in accordance with the fol- 
lowing physical picture. The term dudxH appears in 
the energy current expression thereby whenever A = 
the massless Higgs field do transport energy to 
The extent of this energy transport is restricted only 
by the initial data and the field equation. 

Interesting subcases occur (which would deserve 
further investigations) whenever 

• A = and i/^ = 

• A 7^ and i/o 7^ but iJ„ = 

• A 7^ and i/o = 0. 

Then either or both of the fields are massless and 
there is no restriction on the behavior of the actual 
massless field or fields at J^+. 

The above considerations shows that a YMH sys- 
tem satisfying reasonable regularity assumptions at 
cannot radiate energy to future null infinity un- 
less either A = 0, i?o = or =0. 



Appendix B 

This appendix is to list expressions we applied to 
represent various numerical i?-derivatives in miscella- 
neous stencils. The expressions below always refer to 
an arbitrary function / on the Z*'' time slice. 



31 



The first order i?-derivative relevant for a symmet- 
ric fourth order stencil reads as 



1 



12AR 



( fi+2 



- 8/, 



i-l 



fi-2) 



(B.l) 

The sixth order i?-derivative relevant for a symmet- 
ric sixth order stencil was approximated as 



(B.2) 



We also used a numerical adaptation what would be 
called as 'one sided derivatives' in analytic investiga- 
tions. The relevant first order 'left sided derivatives' 
are approximated in our fourth order stencil as 



-16/1 



1 



12AR 



(3/1 



(B.3) 



48/j 



25/1, 



and as 



+6/1 



1 



\2AR 

-3 - 18/L.. 



-/I 



(B.4) 



io/l„ 



3/L.. 



The analogous sixth order 'left sided' derivatives in 
a sixth order stencil at the spatial gridpoints with in- 
dices Imax-, I-max - 1 and I^ax - 3 do uot differ from 
each other. They are simply yielded by the substitu- 
tion oil = Imax — 3 into i|B.2p . 
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